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A  study  is  made  of  two  basically  distinct  methods  normally  used 
in  the  development  of  band  models  for  the  calculation  of  gaseous 
olecular  transmittance  in  the  infrared  region.  The  first  method 
onsists  of  the  determination  of  the  fempirica  ir  transmittance 
unction  and  the  associated  absorber  and  spectral  parameters  from 
easured  or  calculated  transmittance  spectra.  The  second  method 
onsists  of  the  determination  of  the  absorber  and  spectral  parameter! 
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with  an  assumed  ^analytical*  transmittance  function,  using 
the  same  type  of  data.  Computerized  numerical  techniques 
are  presented  in  connection  with  the  first  method  and  a 
generalized  transmittance  function  is  adopted  for  the  second 
method.  Although  the  methodology  is  generally  applicable 
to  other  gaseous  species,  it  is  specifically  discussed  in 
connection  with  the  trace  gases  SOjf,  NO,  NO?,  and  NHtf.  As 
a  secondary  '  ef f or t  a  structural  breakdown  or  the  Lowtran 
code  is  presented  for  the  purpose  of  incorporating  the  band 
models  for  the  trace  gases.  The  code  is  separated  into  basic 
functional  modules  or  subroutines  controlled  by  a  main  pro¬ 
gram.  The  modularization  itself  was  primarily  performed 
under  a  separate  effort  through  the  Atmospheric  Sciences 
Laboratory.  \ 


As  eM6ioaV»r 

ms  mm 

xdc  rs 

ttMMMHUetd 

Justification. 


^lstrlbu*  Vrn / . — 

AvritnblLlt7./'  ■,,5. 
kvaiflB'd/oi- 


TABLE  OF  CONTENTS 


Page 

Abstract . ,  lv 

List  of  Tables, . . . . . .  vill 

List  of  Figures . . . .  . .  x 

I.  Introduction.., . 1 

II.  The  Transmittance.  Equation .  3 

III.  Computerized  Method  of  Empirical  Model 

Development .  . .  6 

3.1  Introduction, . 6 

3.2  Data  Structure...., . 7 

3.3  Linearization  of  the  Transmission 

Functions . 7 

3.4  Formation  of  the  Square  Error.,,..,...,  9 

3.5  Auxiliary  Variables...., . 10 

3.6  Regression  Analysis.., .  11 

3.7  Piecewise  Analytical  Transmittance 

Function, . . . 14 

e 

H 

3.8  Cf  for  Non-major  Bands,, . .  ......  16 

IV.  Computerized  Method  of  Analytical  Model 

Development . 17 

4.1  Introduction .  17 

4.2  Basic  Equations. . 17 

4.3  Nonlinear  Optimization  Method .  19 

V.  Comparison  of  the  Two  Methods.,,,,,, .  22 

5.1  Introduction., .  22 

5.2  Final  Products, . 22 

5.3  Installation  of  Results  in  Lowtran.  .  .  .  .  23 


ill 


Table  of  Contents 


(continued) 


Page 


5.4  Data  Requirements . 23 

5.5  Computation  Time.,...,.,,,,..,..,.......  24 

VI.  Lowtran  Capabilities  and  Functions .  27 

6.1  Introduction... . 27 

6.2  General  Features, . 29 

6.3  Resonant  Molecular  Absorption  Models....  35 

6.4  Non-Resonant  Molecular  Absorption 

Models  37 

6.5  Scattering  Models. .  42 

VII.  Modularization  of  Lowtran  Including  the 

Trace  Gases . 43 

7.1  Introduction . 43 

7.2  Structure  of  Modularized  Version .  43 

7.3  Models  for  H20  Vapor,  Infrared  03 

and  the  Unif ormly-Mlzed  Gases .  55 

7.4  Models  for  Trace  Gases  SO-,  NO,  N02, 

and  NH3  .......  .  56 

VIII.  Calculations  and  Results,...,, .  75 

8.1  Introduction...,, . 75 

8.2  Testing  of  Modularized  Version .  76 

8.3  Band  Model  Development . . 86 

IX.  Discussion  and  Conclusions . . .  93 

9.1  Introduction.,.,.,...,,.,', .  93 

9 . 2  Changes  . . . .  •  •  94 

9.3  Model  Development,.,,,. . . . 95 


iv 


Table  of  Contents  (continued) 

Page 

9,4  Data  Structure...,,.,..,..., .  98 

References . . . . . .  101 

Appendix:  Modularized  Computer  Code  with  the 

Trace  Gases . . .  103 

Programs:  ADSET .  155 


SIMMIN 


173 


LIST  OF  TABLES 


Page 


Frequency  interval  of  the  attenuation  band 
models  in  the  Lovtran  code.  The  models  for 
aerosol  absorption  and  aerosol  and  molecular 
scattering  are  spectrally'-continuous  and, 
therefore,  not  shown . . 


Spectral  region  over  which  the  attenuation 
models  in  Lowtran  are  effective, . 


Absorber  parameters  in  Lowtran  for  the 
attenuation  models,  where  p  is  the  density, 
Z  the  range  and  M  the  mixing  ratio.  The 
H~0  continuum  model  is  excluded  because  of 
its  different  functional  form . 


Concentration  of  the  unif ormly-mixed  gases 
used  in  the  combined  model . 


Constants  for  the  curve  segments  in  the 
empirical  transmission  function  for  S02 , 

Constants  for  the  curve  segments  in  the 
empirical  transmission  function  for  NO. . 


Constants  for  the  curve  segments  in  the 
empirical  transmission  function  for  NO^. 

Constants  for  the  curve  segments  in  the 
empirical  transmission  function  for  NH^ . 


Absorber  parameters  in  Modularized  Lowtran 
used  with  the  models  for  the  trace  gases... 


The  spectral  coefficient  C’(v)  for  S02. 
The  spectral  coefficient  C'(v)  for  NO.. 
The  spectral  coefficient  C'(v)  for  N02. 
The  spectral  coefficient  C'(v)  for  NH^ , 


Absorption  frequency  region  of  the  trace 
gases  in  the  atmosphere . 


vi 


List  of  Tables  (continued) 


?age 

9.  Transmittance  difference  between  calculations 

using  the  tabulation  of  the  transmittance 
functions  and  calculations  using  the  continuous 
function  representation  for  a  5  km  path  at 
sea  level  in  a  subarctic  winter  atmosphere....  77 

10.  Calculations  of  trace  gas  (T.G.)  trans¬ 

mittance  for  a  5  km  path  at  sea  level  In  a 
tropical  atmosphere  with  a  23  km  visual 
range.  The  columns  on  total  transmittance 
include  all  the  attenuations  and  the  trace 
gases,  except  for  the  rightmost  column  which 


excludes  the  trace  gases .  80 

11.  a.  Atmospheric  regions  Included  in  the  data 

calculations, . 81 

b.  Transmittance  cuts  chosen  from  the  curve 

of  growth .  81 

12.  a.  Band  model  parameters  for  SO^ . 82 

b.  Band  model  parameters  for  NO . 83 

c.  Band  model  parameters  for  NOj . . .  84 

d.  Band  model  parameters  for  NH^  .  .  .  .  . .  85 


13,  Comparison  of  standard  deviations  in  t  .  The 

two  rows  ADSET  1  and  2  are,  respectively,  for 
the  piecewise  analytical  transmission 
functions  with  linear  and  quadratic 
exponents .  92 


LIST  OF  FIGURES 


Page 


Schematic  representation  o£  "equal 
transmittance"  data  structure...... 


Schematic  representation  of  linear  shift 
accounting  for  spectral  dependence  of 
transmittance . .  .  - . . 


Pre-processing  of  data  (t.,  x.),  and  the 
derived  equal  transmittance  data . 


Conceptual  flow  chart  of  modularized  Lovtran. . 

General  flow  chart  for  modularized  Lowtran  4.. 

Flow  chart  for  subroutine  ABSORB,  computing 
equivalent  path  quantities . 


Flow  chart  for  transmittance/radiance  loop.... 

Flow  chart  for  subroutine  FREQSL  and 
thirteen  transmittance  model  subroutines . 


9.  a.  Empirical  transmission  function  for  SO2 


Empirical  transmission  function  for  NO. 


Empirical  transmission  function  for  N02* 
Empirical  transmission  function  for  NH^. 


10.  a.  Standard  transmission  function  from  ADSET 


with  a^=0, 


b.  Standard  transmission  function  from  ADSET 


with  a^O. 


Standard  transmission  function  from  SIMMIN... 


itmmi 


1 


I.  Introduction 

Following  the  efforts  of  Elssaser*  numerous  workers 
have  attempted  to  arrive  at  computationally-simple  models 
for  gaseous  molecular  transmittance,  averaged  over  narrow 
spectral  intervals  in  the  infrared.  These  efforts  may  be 
naturally  divided  into  those  involving  the  analytical 
derivation  of  a  mean  transmittance  function  from  Beer’s 
Law,  and  those  involving  the  extraction  of  the  transmittance 
function  itself  from  transmittance  data.  Traditionally, 
the  former  are  called  "analytical"  and  the  latter  are  called 
"empirical".  The  method  normally  used  in  the  empirical 
models  consists  of  the  extraction  of  the  transmittance 
function  through  graphical  techniques,  with  the  adoption 
of  a  relationship  between  spectral  and  absorber  parameters. 
In  the  development  of  analytical  models  a  transmittance 
function  is  adopted  at  the  offset,  and  the  spectral  and 
absorber  parameters  are  afterward  determined  through  com¬ 
puterized  numerical  procedures. 

In  the  work  reported  here  the  authors  present  a 
totally  computerized  version  of  the  classical  graphical 
methods  for  the  extraction  of  the  empirical  transmittance 
function.  This  is  followed  by  a  presentation  of  a  numerical 
method  which  uses  a  double-exponential  transmittance  func¬ 
tion  for  the  development  of  analytical  band  models.  Both 
methods  are  then  applied  to  20  cm 


averaged  line-by-line 


2 


transmittance  data  for  the  atmospheric  trace  gases  SC^* 

NO,  NC^,  and  NH^.  The  model  parameters  are  listed  at 
5  cm  ^  intervals  throughout  the  major  absorption  bands  of 
these  gases  for  the  convenience  of  the  community  of  band 
model  users.  Although  the  methodology  is  applied  speci¬ 
fically  to  the  trace  gases,  no  restrictions  are  immediately 
evident  in  the  extension  to  other  gaseous  absorbers  in  the 

infrared.  In  fact,  the  analytical  method  was  successfully 
2 

applied  earlier  to  the  principal  band  centers  of  the  major 
absorbers  1^0  vapor,  0^  and  the  uniformly-mixed  gases. 

As  an  application  of  the  results  found  through  this 
effort,  the  band  models  for  the  trace  gases  were  incor¬ 
porated  in  the  widely-used  code  called  Lowtran.  To 
facilitate  the  inclusion  of  these  models,  as  well  as  of 
others,  the  code  was  broken  down  into  separate  subroutines 
or  modules  controlled  by  a  master  program.  The  subroutines 
include  the  evaluation  of  the  equivalent  absorber  amount, 
the  selection  of  the  spectrally-ef f ective  attenuation  model 
and  the  individual  attenuation  models.  The  principal 
purpose  of  the  modularization  is  to  assist  users  with  the 
modification  of  the  code  to  suit  their  individual  re- 


> 


quirements  on  transmission  models. 


/ 
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II.  The  Transmittance  Equation 

The  monochromatic  transmittance  at  frequency  V 
for  the  passage  of  infrared  radiation  through  a  path  length 
Z  in  an  inhomogeneous  medium  with  pressure  and  temperature 
distributions  P(Z)  and  T(Z),  respectively,  is  given  by 
Beer's  Law  in  the  form 

%  .  "I  VP-T>d"<Z>  (1) 

where  is  the  resultant  absorption  coefficient  for  all 
contributing  lines  and  gaseous  absorbers,  and  U  is  the 
absorber  amount.  For  broadband  radiation  detected  by  an 
instrument  of  spectral  response  the  variable  of  interest 

is  the  weighted  mean  transmittance  T,  defined  as 


T 


J  Vvdv  ' 


|  V 


(2) 


Equation  (2)  has  been  evaluated  analytically  over  a 

spectral  interval  Av  for  the  special  case  of  Lorentzian 

broadened  lines  having  assumed  line  distributions  and 

1  3 

intensities,  leading  to  the  classical  band  models  ’  . 

Numerous  variations  of  the  classical  band  models  may  be 
found  in  the  literature,  most  of  which  specify  the  analytical 
form  of  t  in  terms  of  mean  line  or  meteorological  variables. 

4 

A  notable  exception  is  the  model  of  King  which  expresses 


the  homogeneous-path  transmittance  as 


g (Sa  U) , 


(3) 


where  g  is  a  function  to  be  determined  empirically,  S  is 
the  mean  line  intensity,  a  is  the  mean  line  half-width  and 
n  is  an  absorber  parameter  with  the  physical  constraints 
of  zero  and  one  in  the  weak-line  and  strong-line  limits, 
respectively.  The  path  inhomogeneity  may  be  accounted  for 
in  Eq .  (3)  through  the  Cur t is-God son  equivalences 


SanU  =  j 


S(Z)an(Z)dU(Z)  . 


(4) 


From  practical  considerations,  it  is  often  desirable 
to  transform  the  argument  in  Eq.(3)  with  the  known  relations 


S 


a 


T  a 
So  <T*> 


T  h 

„  w  o.1 

o 


(5) 

(6) 


in  order  to  obtain 


T  *  g  {  C(f-)  c^)  U  }  ,  (7) 

o 

where  C  is  -a  spectral  parameter  combining  Sq  and  a”,  m  is 
an  absorber  parameter  combining  the  temperature  exponents 
of  S  and  a,  a  is  an  absorber  constant,  and  the  subscript 
"o"  denotes  standard  conditions.  For  computational  con¬ 


venience  Eq.  (7)  may  be  expressed  as 


5 


where 

x  -  C»  +  log1QW 


(9) 


C*  -  log1QC  (10) 

n  T  m 

W  -  Cf— )  (^)  U  .  (11) 

1  o 

Here,  £  Is  the  transmittance  function,  C'  Is  the  spectral 
parameter,  W  Is  the  equivalent  absorber  amount,  and  n  and 
m  are  the  absorber  parameters;  all  of  which  are  to  be 
determined  from  transmittance  data  for  each  absorber. 
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III.  Computerized  Method  of  Empirical  Model  Development 
3.1  Introduction 

Assuming  the  availability  of  equal  transmittance 
data,  which  is  defined  below,  we  have  developed  an  algorithm, 
called  ADSET,  which  evaluates  absorber  parameters  n,  m, 
spectral  parameters  C ’  (v )  and  an  empirical  transmission 
function  simultaneously.  In  the  algorithm  the  transmission 
function  is  linearized  and  a  linear  regression  technique 
is  utilized  for  parameter  evaluation.  In  order  to  evaluate 
the  band  model  parameters  and  the  empirical  transmission 
function  simultaneously,  a  set  of  auxiliary  variables  are 
introduced.  Each  data  point  is  identified  through  the 
auxiliary  variables  to  an  absorption  band  and  to  a  trans¬ 
mittance  'cut'.  This  enables  us  to  obtain  globally  optimal 
set  of  parameters  and  the  empirical  transmission  function 
simultaneously . 

Based  on  the  derived  optimal  pointwise  transmission 
function,  a  piecewise  analytical  transmission  function  is 
developed.  The  commonly  used  computer  code  Lowtran  for 
the  evaluation  of  atmospheric  transmittance  can  be  greatly 
simplified  by  the  use  of  this  piecewise  analytical  trans¬ 
mission  function  to  model  the  major  absorbers. 

Finally,  the  code  ADSET  also  contains  a  subroutine 
which  can  compute  the  spectral  parameter  value  C'(v)  for 
non-major  absorption  bands. 


3.2 


Data  Structure 


Several  transmittance  values  Tj ,  j«l,  2,  .  ..,  NCUT 
are  chosen  a  priori,  where  NCUT  is  the  number  of  chosen 
transmittance  values.  Curves  of  growth  data  (l.e.  T  versus 
U)  for  each  layer  of  atmosphere  are  assumed  to  be  given  at 
these  transmittance  values.  Therefore,  the  curves  of 
growth  have  'cut'  structure,  namely,  all  data  points  are 
on  one  of  the  transmittance  cuts  T  ■  Tj  ,  j“l,  2,  ...,  NCUT 
(.See  Fig.  1).  We  call  a  data  set  with  this  cut  structure 
an  'equal  transmittance'  data  in  the  sequel. 

3.3  Linearization  of  Transmission  Function 

Since  f  in  Eq.(8)  is  known  to  be  strictly  monotone 
decreasing  from  one  to  zero  as  x  changes  from  to 
there  exists  an  inverse  function  f  ^  defined  on  (0,1)  such 
that 

x  *  f  ^  (t  ) 

-  c'  +  logW 

T 

-  C'  +  nlog(|-)  +  mlog^)  +  logU.  (12) 

o 

Let  us  define  x^  ,  j=*l,  2,  ...,  NCUT  be  the  inverse  image 
of  the  prechosen  transmittance  values  Tj ,  j*l,  2,  ..., 

NCUT  i.e., 

X5  "  f  1CV’  J“  1>  2 . 


NCUT, 


(13) 


Absorber  Aiaount 

(logarithmic  scaler 
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Then,  the  set  of  points  (xj.Tj),  J  >1,  2,  NCUT  is 

nothing  but  the  empirical  transmission  function,  which  is 
to  be  found.  From  Eq.  (12),  we  reach  the  following  re¬ 
gression  equation. 


n  log (p-)  +  «  logC^)  +  c*  “  x  *  -logU.  (14) 

o 

Note  that  this  equation  is  linear  in  the  unknown  parameters 
n,  m,  C*  and  x.  Therefore,  the  linear  regression  technique 
can  be  used  to  evaluate  the  optimum  values  for  the  para¬ 
meters. 


3.4  Formation  of  the  Square  Error 

The  square  error  corresponding  to  the  k-th  data 
point  in  i-th  absorption  band's  j-th  cut,  denoted  by  E^^, 
is  given  by 


Eijk  *  <«l°8(-f^>  +  mlogOj-M  +  C»  -  Xj  -  (-logUljk)> 


Hence,  the  total  square  error  E^j  for  this  cut  is 


Eij  "  [ml  Eijk,  (16) 

where  is  the  number  of  layers  in  this  cut.  Similarly, 

the  total  square  error  E^  for  1-th  band  and  the  grand 
total  square  error  E  are  given  by 


llj 
E  E , 


j-1  k-1 


NB  NB  Ji  Lij 

E-E  E-E  E  E  E...  ,  (18) 

1-1  1  i-1  j-1  k-1  13 K 

where  and  NB  are  the  numbers  of  the  cuts  in  i-th 
absorption  band  and  of  the  absorption  bands,  respectively. 
The  final  expression  can  be  simplified  if  we  assume  that 
the  number  of  layers  (L^)  in  every  cut  is  equal  to  a 
constant  L^.  For  this  case 

NB  Ji  L± 

E-E  E  E  E  (19) 

i-1  j-1  k= 1  J 


Our  objective  is  to  find  optimum  set  of  parameters  (n*, 

m*.  C^*,  C’2*,  ....  C'NB*.  x^,  x2* . XNCUT*^  which 

minimizes  this  grand  total  error  E. 


3.5  Auxiliary  Variables 


In  order  to  perform  the  minimization  of  E  with 
respect  to  the  above  parameters  simultaneously,  we  modify 
the  square  error  E^^  80  that  it  contains  all  the  parameters 
This  is  done  by  Introducing  two  sets  of  auxiliary  variables 

u^,  i-1,  2 . NB  and  v^ ,  j-1,  2,  ...»  NCUT.  Using  them, 

E^k  Is  redefined  as 


!ijk  "  {nl°8(  Pq  >  +  ml°8^T1°k)  +  U1  ,'i  j k  C1  +  +  uNB,ijkC  1 

2 

+  vl,ijkKl  +  +  VNCUT,ijk  ScUT  r  ^"lo8Uijk^  (20) 
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where  ■  -x^ ,  j“l#  2»  .•.,  NCUT.  The  auxiliary  variables  act  as 
identifiers  of  the  band  and  the  cut.  If  a  data  point  is  for  1-th 
band's  1-th  cut,  then  u^  -  1  and  u^  -  0  for  all  i  4  1  and  -  1 

and  Vj  ■  0  for  all  j  4  }.  Thus,  only  the  spectral  parameter  C'^ 
and  the  cut  parameter  Kj  corresponding  to  the  current  data  are  active 
and  all  other  spectral  and  cut  parameters  disappear.  Hence,  Eq.  (20) 
reduces  to  Eq.  (15).  The  change  from  to  ■  -x^  is  made  in  order 
to  symmetrize  the  coefficient  matrix  of  the  resulting  normal  equation. 
This  change  makes  it  possible  to  utilize  any  specialized  solution 
method  for  the  symmetric  normal  equation  when  the  space  conservation 
is  important. 

3.6  Regression  Analysis 

Using  the  grand  total  error  E  with  the  redefined  in 

It  It  &  ft 

Eq,  (20),  the  best  parameter  values  n  ,  m  ,  C'^  ,  ...,  C'NB  , -K^,..., 

it 

KNCUT  are  3ilnultaneoualy  determined  by  the  linear  regression.  Setting 
the  partial  derivatives  of  E  with  respect  to  parameters  equal  to 
zero  results  in  a  linear  normal  equation  of  the  form  AX  “  B,  where 
A,  B  and  X  are,  respectively,  a  symmetric  coefficient  matrix,  a 
constant  vector  and  a  parameter  vector  defined  by 


Eu2  n 

NB  0 

1  A  “ 

* 

• 

* 

0  * • Su22 
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B 


jl(-vNCUT  logU),  ....  EC^v^ogU),  Z(-uNglogU)  ,  .  .., 

T  t 

Z(-u2logU),  I(-log(~)logU),  S(«-log(|-)logu2  ,  (22) 


£ 


N  C  U  T  ' 


r  » 

CNB’ 


C  1 
C2 


m,nj 


(23) 


The  *  in  Eq .  (21)  represents  some  nonzero  elements.  Also 

the  £  in  the  above  equations  represents  the  triple  sum 
NB  J±  L 

£  £  £  in  Eq.  (19).  One  may  realize  that  C'  does  not 

i=l  j=l  k-1 

appear  in  Eq .  (23)  and  hence  the  corresponding  auxiliary 

variable  u^  is  also  absent  from  Eqs.  (21)  and  (22).  This 

is  because  one  of  C’ ,  ....  C '  is  dependent  on  other  C ' 

1  N  J5  1 

so  that  C*  ....  C *  cannot  be  determined  uniquely.  It  is 

1  N  D 

necessary  that  one  of  C^s  be  given  a  number  a  priori.  Here 
is  chosen  and  is  given  the  value  zero,  and  therefore  is 
eliminated  from  the  parameter  vector  X.  This  choice  calls 
for  some  explanation.  On  t  vs.  logW  diagram  the  optimum 
empirical  transmission  function  can  be  placed  anywhere. 

What  it  amounts  to  is  that  a  different  placement  results 
in  a  different  set  of  values  which  is  a  linear  shift 
(addition  or  subtraction  of  a  constant)  of  another  set  of 
C|  values.  Only  the  relative  relationship  among  C|  is 
unique.  This  is  clearly  indicated  in  Fig.  2. 

Since  the  placement  of  the  empirical  transmission 


function  is  arbitrary,  we  may  position  it  on  the  data 
points  corresponding  to  the  first  absorption  band.  In 
other  words,  the  first  absorption  band  is  taken  as  the 
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reference  band  and  the  corresponding  spectral  parameter 
Cj  is  set  to  be  zero. 

The  queuing  of  parameters  in  X  vector  is  determined  in 
such  a  way  that  as  many  upper  principal  minor  matrices  as 
possible  become  diagonal  (See  Eq.  (21)).  This  arrangement 
can  reduce  the  amount  of  computation  in  the  early  stage  of 
Gauss  elimination  steps  when  the  normal  equation  is  solved, 
and  can  result  in  less  computational  error. 


3.7  Piecewise  Analytical  Transmission  Function 

After  the  best  parameter  values  are  computed,  the 
piecewise  analytical  transmission  function  is  generated 
by  the  piecewise  interpolation.  The  transmittance  region 
(0,1)  is  divided  into  NCUT  -  1  subregions  by  the  trans¬ 
mittance  cuts  x2>  T  ^  >  •••,  tncUT-1‘  Let  T1  >  X2  >  •••  > 

TNCUT’  t*ien  the  subregions  are  given  by 


Subregion  1 
Subregion  2 


[t2,1  ), 
[x3 ,1^]  , 


Subregion  NCUT-2 
Subregion  NCUT-1 


^TNCUT-l,TNCUT-2^  ’ 
^ ' TNCUT-1 ^ * 


The  top  and  bottom  subregions  contain  T ^  and  t 


NCUT 


as  an  inner  point,  respectively.  The  interpolation  in  each 


subregion  is  done  by  the  double  exponential  function  defined 
by 

2 

a.  +  a_x  +  a.x 

t(x)  =  exp  {-10  1  )  .  (24) 

The  generally-used  linear  interpolation  is  not  used  here 
since  subregions  cannot  be  assumed  small  enough  for  the 
linear  approximation  to  be  valid.  Furthermore,  the  linear 
interpolation  is  totally  inadequate  for  the  top  and  bottom 
subregions.  On  the  other  hand,  the  double  exponential 
function  takes  the  values  between  and  is  asymptotic  to  one 
and  zero  as  the  argument  varies  from  to  00 .  It  is  also 
known  that  this  function  closely  approximates  the  standard 

2  i 

empirical  transmission  function  used  in  the  Lowtran  code.  ’ 
The  parameters  a^ ,  and  a^  for  each  subregion 

are  determined  by  two  different  methods.  The  first  method 
assumes  that  a^  =  0  and  uses  no  further  data  to  compute 
a^  and  a They  are  simply  determined  by  the  condition  that 
the  interpolation  function  in  each  subregion  passes  through 
the  end  points.  In  the  top  and  bottom  subregions,  the  func¬ 
tion  is  required  to  pass  through  two  points;  (t^,  x1^  anc* 

(t2,  x 2 )  for  the  top  and  (tnCut-1  ,XNCUT-1^  and  ^NCUT’ 
for  the  bottom  subregions. 

NOU  i 

The  second  method  does  not  assume  that  a^  *  0  and 
requires  additional  data  to  compute  parameter  values.  The 
same  condition  that  each  Interpolation  function  passes 


through  two  points  reduces  the  number  of  unknown  parameters 
to  one.  The  last  parameter  is  determined  by  minimizing  the 


subregional  square  error  defined  by 

L.  .  2 

i  a.  =  a0x.  +  a,x. 

E  =  Z  (T  -  exp  {-10  1  3  )) 

i=l 


(25) 


for  those  data  points  in  respective  subregions. 

3.8  C'  for  Non-major  Bands 

Finally,  the  spectral  parameters  C*  for  non-major 
bands  are  computed  by  a  straightforward  method.  The  dis¬ 
crepancies  between  x*j,  and  logV^  values  computed  for  all 
cuts  for  one  band  are  averaged  to  obtain  the  spectral 
parameter  C'(v)  for  that  band,  i.e., 

1  N 

C'  =  i  E  (x*  -  logW  ) ,  (26) 

i  =  l 

where  are  computed  by  Eq .  (11)  with  optimal  n*  and  m* . 


IV. 


Computerized  Method  of  Analytical  Model  Development 


4.1  Introduction 

In  the  last  chapter,  we  assumed  no  analytical  form 
for  the  transmission  function  T  =  f  (x)  when  the  standard 
transmission  function  was  computed.  Here,  by  assuming  the 
double  exponential  form  given  by  Eq.  (24)  as  the  transmission 
function  for  the  entire  transmittance  range,  we  derive  an 
algorithm  which  can  evaluate  the  best  function  parameter 
values  a^,  a 2,  and  a^;  together  with  the  band  model  para¬ 
meters  n,  m,  and  C^.  Note  that  the  double  exponential 
function  was  used  for  the  piecewise  interpolation  in  the 
last  chapter.  But  the  computation  of  the  function  para¬ 
meters  was  performed  after  the  band  model  parameters  and 
the  empirical  transmission  function  were  obtained.  In 
other  words,  the  computation  in  the  last  chapter  was 
sequential  but  not  simultaneous.  The  algorithm  we  present 
in  this  chapter  is,  on  the  contrary,  the  simultaneous 
evaluation  of  all  parameters.  The  preliminary  development 
of  this  algorithm  can  be  found  in  Ref.  5. 

4.2  Basic  Equations 

The  basic  equations  are  Eq.  (8)  and  Eq.  (24)  of 
the  last  chapter,  which  are  cited  here  for  the  ease  of 
reference . 

T  -  f(x)  2  (27) 

a  +  a_x  +  a,x 

f(x)  -  exp  {-10  1  }.  (28) 


jfoatiw&stff  Tipr 


Now,  since  we  have  assumed  the  function  form,  we 
can  compute  the  transmittance  if  we  have  the  value  of  x. 
Hence,  we  do  not  have  to  take  the  inverse  function  as  we 
did  before  to  perform  the  regression  analysis.  Instead, 
we  take  the  square  difference  of  the  given  and  computed  T 
directly  from  this  expression.  Thus,  we  get 

2  2 

8,+a.x  .+a_x . . 

E  =  [T±j  -  exp { -10  1  3  3  }]  ,  (29) 

for  i-th  absorption  band's  j-th  data  point,  where,  as 
before,  x^  is  given  by 


r  .  .  T 

=  C|  +  n  log  ( p~^~ )  +  logfcp-)  +  logU  . 


By  summing  this  individual  error  for  all  data  in  i-th  band, 
we  have  the  total  error  for  this  band  as 


E  =  E  E  (31) 

j  =  l  J 

where  is  the  number  of  data  in  i-th  band. 

Again,  we  introduce  auxiliary  variables  u^,  i=l,2, 

. . . ,NB  in  order  to  introduce  all  C^,  i=l,2,...,NB  into  the 
x^j  expression  Eq.(30).  By  this  we  get 

NB  P  T 

xij  =  l  ,  Uk,ijC'k  +  "  l08(P  >  +  m  l08(T77)  +  l0gUij’  ( 


We  use  this  expression  for  x . .  in  the  following  total  error  E 


l 


ii  •ijMwfr'  r*  « 


NB 

E  =  Z 


i  =  l  j  =  1 


Now,  we  are  ready  to  take  partial  derivatives  with 
respect  to  the  parameters  n,  m,  C ' ,  C '  ,  a, ,  a _ ,  and 

a^  to  form  the  normal  equation  for  this  regression  problem. 
Theoretically  speaking,  we  can  evaluate  the  'best*  para¬ 
meter  values  by  solving  the  normal  equation.  But  obviously 
the  grand  total  Eq .  (33),  which  is  to  be  minimized,  is  not  a 
quadratic  function  of  the  unknown  parameters  and,  therefore, 
the  resulting  normal  equation  is  not  a  linear  function  of 
them.  Hence,  we  need  to  adopt  a  different  numerical 
method  for  the  evaluation  of  the  'optimal'  parameter  values. 


4.3  Nonlinear  Optimization  Method 


Th.e  computational  technique  we  use  here  is  a  recursive 
technique  which  is  referred  to  as  the  conjugate  gradient 
method^.  jn  essence,  this  technique  improves  a  set  of 
guesses  of  the  parameter  values  recursively  by  locating  a 
new  set  of  guesses  which  yields  smaller  error.  For  a  given 
guess  (an,3n,  ...,yn)  of  the  minimizing  parameter  vector, 

at  which  the  error  is  minimized,  the  best  direction  of  the 
search  in  the  parameter  space  for  a  new  guess  is  first 
determined  using  up  to  second  order  derivatives  of  the 
error.  Then  the  one-dimensional  search  for  the  minimizing 
point  is  performed  along  this  direction  from  (an,gn,  .... 
y  )  to  find  a  new  guess  (ex  ,p  ,  ...,Y  ) 
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+  .2x}  +  a3xj 


£U1>' 


and  f  (x)  is  given  by  Eq.  (28). 

Note  that  there  exists  a  linear  dependence  among  the 
gradients  which  is 


Therefore,  the  parameter  set  {n,  m, 
a^)  cannot  be  determined  uniquely. 


C1 . CNB  ’  al»  a2  * 

As  it  was  explained  in 


the  previous  section,  this  is  due  to  the  arbitrariness  in 


the  positioning  of  the  standard  transmission  function. 


Hence,  the  spectral  parameter  C|  is  again  set  to  be  zero, 
so  that  we  can  evaluate  unique  set  of  optimal  parameters. 
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V.  Comparison  of  the  Two  Methods 

5.1  Introduction 

Both  methods  can  evaluate  the  optimal  n,  m  and  C'(v) 
values  for  major  and  non-major  bands  and  also  a  standard 
transmission  function.  But  there  are  some  basic  differences 
which  are  discussed  in  the  sequel. 

5.2  Final  Products 

The  final  product  of  the  ADSET  code  is  a  piecewise 

analytical  standard  transmission  function  together  with 
the  band  model  parameters.  Each  analytical  piece  of  the 
standard  transmission  function  covers  only  one  of  the  pre¬ 
chosen  subintervals  of  (0,1)  transmittance  range.  On  the 
other  hand,  SIMMIN  produces  only  one  analytical  trans¬ 
mission  curve  for  the  entire  range.  Therefore,  ADSET  has 
more  flexibility  to  adjust  to  the  transmittance  curve 
variations.  This  feature  of  ADSET  can  be  very  valuable 
for  the  gases  with  non-standard  curves  of  growth. 

This  difference  is  amplified  when  the  number  of  the 
transmittance  sub-regions  used  in  ADSET  is  increased. 
However,  as  the  number  of  subregions  increase,  the  re¬ 
quirement  on  the  usable  data  becomes  severer  and  more 
spaces  are  necessary  to  store  the  computed  results.  Hence, 
the  determination  of  the  number  of  subregions  should  be 


resorted  to  compromise. 


5.3 


Installation  of  the  Results  in  Lowtran 


The  final  products  of  two  codes  ADSET  and  SIMMIN 
were  installed  in  the  widely-used  Lowtran  code,  as  discussed 
in  Section  VII  of  this  report.  The  SIMMIN  results  require 
less  memory  space,  less  time  for  transmittance  computation 
and  simpler  coding  than  the  ADSET  result.  In  fact,  for  the 
SIMMIN  result,  all  that  have  to  be  stored  are  the  five 
band  model  parameters  n,  m,  a^ ,  a2»  and  a^,  and  a  set  of 
spectral  parameters  C'(v^)  for  each  absorber.  Furthermore, 
the  computation  of  t  can  be  done  by  only  one  FORTRAN 
statement.  On  the  other  hand,  the  ADSET  result  requires 
the  storage  of  NCUT-1  of  a^  ,  a^t  and  values,  n  and  m 
and  a  set  of  C'(v)  values  for  each  absorber.  There  can  be 
a  large  difference  in  the  number  of  the  sets  of  a  ,  a.}, 
and  a  ^  values  to  be  stored.  Moreover,  some  judging  state¬ 
ments  are  necessary  to  select  the  right  set  of  a^,  a,,  and 
a^  for  each  transmittance  computation. 

5.4  Data  Requirements 

The' ADSET  code  requires  the  cut  structured  data 
such  that  the  transmittance  of  each  data  point  must  .  i' 
in  one  of  prechosen  values.  But  the  SIMMIN  code  does  not 
impose  any  conditions  on  the  data  set. 

Some  considerations  on  the  requirement  of  equal 
transmittance  data  for  ADSET  are  due  here.  Even  if  the 
available  data  do  not  have  equal  transmittance  structure. 
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it  can  be  transformed  into  the  required  form  using  inter¬ 
polation/extrapolation.  This  constitutes  the  pre-processing 
of  the  raw  data.  Curves  of  growth  data  x  vs.  log  U  with 
x  values  not  necessarily  coinciding  with  the  prechosen 
values  can  be  locally  interpolated/extrapolated  using  an 
analytical  function.  This  procedure  is  indicated  in  Fig. 

3.  Again,  the  double  exponential  function  is  an  excellent 
choice  for  the  interpolation  function.  We  note  that  an 
almost  exact  technique  as  the  one  ustd  in  obtaining  a 
piecewise  analytical  transmission  function  can  be  used  for 
this  purpose.  In  fact,  only  a  minor  modification  of  the 
interpolation  subroutine  used  in  ADSET  can  accomplish  this 
task. 

5.5  Computation  Time 

The  numerical  methods  used  in  ADSET  and  SIMMIN  for 
solving  normal  equations  are  essentially  different.  The 
method  in  SIMMIN  is  a  recursive  algorithm  and  the  other  in 
ADSET  is  a  non-iterative  one.  Therefore,  the  computation 
time  for  ADSET  is  determined  by  the  size  of  the  data  set 
only,  whereas,  the  one  for  SIMMIN  depends  on  both  the 
actual  data  values  and  the  initial  guesses.  It  is  difficult 
to  estimate  the  computation  time  for  SIMMIN  due  to  this 
dependence.  One  way  of  controlling  the  time  is  to  limit 
the  number  of  iterations  performed.  This  feature  is  included 
in  the  packaged  subroutine  FMCG  which  is  used  for  actual 


VI.  Lowtran  Capabilities  and  Functions 

6.1  Introduction 

The  Lowtran  code  consists  of  a  computer  model  for 
the  calculation  of  transmittance  through  atmospheres 
containing  absorbing  and  scattering  molecules  and  aerosols. 

The  models  used  in  the  code  were  for  the  most  part  developed 
in  1972^  but  later  editions  incorporated  computational 
changes  and  other  capabilities^  It  covers  the  spectral 

range  from  0.25  to  28.5  ym  at  intervals  of  5  cm  *  with  t. 
resolution  for  the  major  absorbers  of  20  cm-1.  The  trans¬ 
mittance  calculation  is  made  on  six  model  atmospheres  and 
two  haze  models  on  a  33-level  basis  for  altitude,  pressure, 
temperature  and  density  from  sea-level  to  100  km.  The  path 
of  the  transmission  is  considered  to  be  refracted  by  changes 
in  atmospheric  density,  a  fact  taken  into  account  in  an 
optional  subroutine.  In  its  present  form  the  Lowtran  code 
consists  of  a  single  main  program  that  inputs  the  path 
data  and  model  parameters,  computes  the  equivalent  absorber 
amount,  and  performs  the  transmittance  calculations.  The 
only  present  subroutines  are  associated  with  the  path,  and 
are  optional.  The  difficulties  of  understanding  and, 
especially,  updating  such  a  program  structure  are  considerable. 

The  principal  objective  of  this  effort  is  to  modify 
the  program  structure  of  the  Lowtran  code  in  accordance 
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with  the  following  criteria: 

1.  The  basic  functions,  calculations  and  print¬ 
outs  remain  nearly  identical  to  the  original. 

2.  The  basic  operations  involving  the  reading  of 
data,  the  calculation  of  the  equivalent  path 
and  the  transmittance  calculations  are  all 
separate,  independent  programs,  but  are  con¬ 
nected  as  subroutines  to  a  main  control  program. 

3.  The  structure  modification  is  performed  on  the 
latest  version  of  the  code,  i.e.,  Lowtran  4. 

As  an  exercise  in  the  use  of  the  modularized  version,  the 
present  authors  added  empirical  band  models  for  trans¬ 
mittance  through  the  trace  gases.  Also,  continuous 
functions  were  made  to  replace  their  transmission  tables 
for  the  principal  molecular  absorbing  species. 
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6.2  General  Features 

In  this  section  an  effort  is  made  to  summarize  the 
basic  structure,  fundamental  calculations  and  models  used 
in  the  Lowtran  code  for  estimating  atmospheric  attenuation 
by  gases  and  aerosols.  Reference  is  specifically  made  to 
the  latest  fourth  version,  although  at  the  present  time 
the  authors  are  aware  of  a  recent  effort  by  AFGL  on  a 
fifth  version.  From  the  authors'  evaluation  of  their 
recent  efforts,  it  appears  that  the  modularization  presented 
here  may  be  incorporated  in  their  latest  version.  For 
instance,  the  latter  is  known  to  have  a  single  separate 
subroutine  for  the  emission  and  transmission  calculations. 
The  major  contribution  of  the  work  presented  here  lies 
in  the  separation  of  that  emission  and  transmission  loop 
into  a  subroutine  for  model  selection,  a  subroutine  for 
the  equivalent  path  and  individual  subroutines  for  all  of 
the  attenuation  models  in  the  code. 

The  Lowtran  code  is  designed  for  the  specific 
purpose  of  calculating  at  low  resolutions  either  atmospheric 
radiance  or  transmittance  between  any  two  locations  in 
the  Earth's  atmosphere  at  frequencies  ranging  from  the 
ultraviolet  (UV)  to  the  infrared  (IR).  This  is  accom¬ 
plished  through  the  use  of  band  models  accounting  for 
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resonant  gaseous  absorption  (e.g.  1^0  vapor,  0^  ,  HNO^ 

vapor  and  the  uniformly-mixed  gases),  resonant  aerosol 
absorption,  non-resonant  gaseous  absorption  (.e.g.  and 

vapor  continua)  and  scattering  by  molecules  and  aerosols 
The  spectral  intervals  over  which  the  band  models  are 
provided  vary  from  5  cm  ^  to  500  cm  ^ ,  as  shown  in  Table  1. 
It  should  be  pointed  out  that  the  spectral  resolution  is 
generally  much  lower  than  the  interval  over  which  they 
are  defined.  For  instance,  the  models  for  the  principal 
absorbers  are  given  at  5  cm  ^  intervals,  while  their- 
spectral  resolution  is  20  cm  1 .  The  spectral  resolutions 
for  the  remaining  models  is  not  specified  anywhere  in  the 
available  literature  on  the  code.  In  this  table  the 
spectral  definition  of  the  models  for  aerosol  absorption, 
and  for  aerosol  and  molecular  scattering  are  not  shown 
because  they  are  spectrally  continuous. 

The  spectral  regions  over  which  the  attenuation 
models  are  effective  are  summarized  in  Table  2.  It  may 
be  seen  In  this  table  that  over  some  regions  only  a  few 
species  attenuate  and,  therefore,  a  transmittance  of  unity 
may  be  specified  in  the  calculation  of  the  total  trans¬ 
mittance.  This  table  forms  the  basis  for  the  model 
selection  subroutine  introduced  in  the  modularized  version 
for  the  purpose  of  simplifying  the  code  structure. 

In  the  discussion  that  follows,  the  individual 
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Table  1.  Frequency  Interval  of  the  attenuation  band  models 
in  the  Lowtran  code.  The  models  for  aerosol 
absorption  and  aerosol  and  molecular  scattering 
are  spectrally  continuous  and,  therefore,  not 
shown . 


ATTENUATING  WAVENUMBER  SCALE  (cm1) 

SPECIES  iOOO  2000  3000  4000  5000  6000  7000 


VIOLET 


35 


attenuation  models  are  grouped  together  in  certain  classes 
and  are  briefly  discussed.  Generally  speaking,  the  dis¬ 
cussion  is  restricted  to  the  extent  of  illustrating  the 
function  and  parameters  which  had  to  be  identified  in 
Lowtran  for  the  modularization  purpose  that  followed.  An 
exception  is  made  in  the  case  of  the  major  molecular 
absorption  models  (i.e.  1^0  vapor,  infrared  0^  and  the 

uniformly-mixed  gases)  because  they  are  replaced  with 
continuous  functions  in  the  modularized  version.  For  a 
comprehensive  discussion  on  the  theory  of  all  of  the 
original  models  the  reader  is  encouraged  to  study  the 
series  of  AFGL  reports^”'1’®  on  the  code,  as  well  as  the 
references  therein. 

6.3  Resonant  Molecular  Absorption  Models 

Molecular  resonant  absorption  is  modeled  in  the 
code  for  b^O  vapor,  infrared  0^,  the  uniformly-mixed  gases 
and  HNO^  vapor.  Different  approaches  are  used  for  the 
first  three  listed  as  compared  with  the  approaches  used 
in  connection  with  0^  in  the  visible  and  ultraviolet 
regions  and  with  HNO^  vapor  in  the  infrared. 

The  models  used  to  account  for  gaseous  absorption 
by  the  molecules  of  HjO  vapor,  infrared  0^»  and  the 
uniformly-mixed  gases  are  based  on  Eq.(8),  namely 

T  -  f  {  x } .  (8  ) 


The  developers  of  Lowtran  obtained  the  parameters  n,  m, 
the  function  f  and  the  spectral  constant  C'  at  5  cm  ^ 
intervals  using  experimental  and  calculated  transmittance 
data  of  20  cm  ^  resolution.  Table  3  shows  the  values  of 
the  parameters,  as  well  as,  the  equations  for  the  calcu¬ 
lation  of  the  absorber  amount.  The  spectral  constant 
C'  over  the  entire  spectrum  of  definition  may  be  found  as 
part  of  the  data  input  presented  in  the  Appendix.  The 
transmission  model  for  the  uniformly-mixed  gases  was 
obtained  by  combining  the  data  for  all  of  these  gases  in 
the  proportions  listed  in  Table  4.  It  should  be  pointed 
out  that  the  temperature  and  pressure  exponents  used  in 
Lowtran  for  the  major  absorbers  and  listed  in  Table  3  are 
not  the  same  as  the  ones  developed  from  the  original 
transmission  data.  This  inconsistency  was  introduced 
during  the  digitizing  of  the  curves  for  inclusion  in  the 

computer  code,  in  order  to  account  more  accurately  for 

2 

the  temperature  dependence  . 

The  method  used  for  modeling  HNO^  vapor  and  the 
visible  and  ultraviolet  0^  is  similar  to  the  one  described 
above  for  the  major  absorbers,  except  that  the  function 
f  was  specified  a  priori  to  be  an  exponential.  Thus, 


T  =  exp  (-CW) , 


where  for  HN0- 
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T 

W  =  (f-)  (~)  U,  (39) 

.  o 

U  =  M  Z  x  105,  (40) 


and  for  0^ 


W  =  U  =  46.667  pZ.  (41) 

In  Eq.  (40)  M  is  the  mixing  ratio  profile  as  tabulated 
in  the  Appendix  together  with  the  C ' s,  and  p  is  the  absorber 
density . 

The  last  of  the  molecular  absorption  models  is  the 
one  for  the  resonant  absorption  by  atmospheric  aerosols. 

The  exponential  function  in  Eq.  (38)  is  assumed 

T  =  exp (-CW) , 


where 

W  =  U  =  3.5336  x  10-6NZ,  (42) 


and  N  is  the  vertical  distribution  of  the  number  of  haze 
particles.  Tabulations  are  provided  of  distributions  for 
5  Km  and  23. km  visibility,  as  listed  in  the  Appendix.  Other 
visibilities  are  treated  in  the  code  itself  through  linear 
interpolation . 

6.4  Non-Resonant  Molecular  Absorption  Models 


Non-resonant  gaseous  molecular  absorption  is  re¬ 
presented  by  the  N„  and  H,0  vapor  continuums.  The  same 
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modeling  approach  is  used  for  ^  as  for  resonant  molecular 
absorption,  that  is 

t  =  exp (-CW) , 

where 

p  2  T  1.5 

W  *  (f-)  U  (43) 

o 

U  =  0.8  Z.  (44) 

For  the  H^O  vapor  continuum  an  exponential  function  is 
also  used,  but  with  a  more  elaborate  exponent.  Thus, 


T  =  exp (-Y) , 


(45) 


where 


C 

Y  -  Cs[Pw  +  ^  (P  -  Pw)]  U.  (46) 

s 

Here,  Pw  is  the  partial  pressure  of  water  and  Cg  and 
are  the  self -broadening  and  nitrogen-broadening  spectral 
constants.  The  values  of  these  spectral  constants  depend 
on  the  spectral  region  where  the  continuum  is  effective. 
In  the  8  to  14  pm  region 

C  =  C  exp  [6.08  -  1)]  ,  (47) 

so  i 


and 


C 

^  =  0.002, 

s 


(48) 


while  in  the  3.5  to  4.2  ym  region 
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C 

s 


C  exp [4.56 
o 


and 


C 

=  0.120. 


In  these  equations  the  value  of  Cq  is  given  by 
Cq  =  4.18  +  5578  exp  (-7 . 87  x  10~3v) . 


(49) 


(50) 


(51) 
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GAS 

MOLECULAR 

WEIGHT 

PARTS  PER  MILLION  BY 
VOLUME  (ppm) 

co2 

44 

330.0 

o 

(Nl 

z 

44 

0.28 

CO 

28 

0.07  5 

CH4 

16 

1.60 

°2 

32 

2.095  x  105 

Table  4 


Concentrations  of  the  uniformly-mixed  gases 
used  In  the  combined  model. 


6.5 


Scattering  Models 


In  order  to  account  for  atmospheric  scattering 
exponential  functions  were  used  again.  For  scattering 
by  molecules  the  model  is  defined  as  in  Eq .  (38) 


T  =  exp(-CW) 


where 


(|~)(^)  U 

(52) 

0 

-2  0 

9.87  x  10  Z 

(53) 

v4 

(54) 

For  aerosol  scattering  the  argument  of  the  assumed  ex¬ 
ponential  function  is 


W  =  U 


3.5336  x  10-4  NZ 


(55) 
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VII.  Modularization  of  Lowtran  Including  the  Trace  Gases 

7.1  Introduction 

Considering  the  generality  and  broadness  in  scope 
of  this  code  it  is  not  surprising  that  the  program  structure 
shows  in  its  present  form  great  complexity.  Although  the 
program  user  is  not  normally  interested  in  aspects  of  the 
code  other  than  the  input  and  output,  there  are  many  cases 
where  a  basic  understanding  helps  in  specific  applications. 
Situations  are  likely  to  occur,  for  instance,  where  a 
replacement  of  one  of  the  several  attenuation  models  is 
highly  desirable.  To  assist  in  the  implementation  of 
model  additions  or  changes  as  well  as  in  the  extension  to 
other  spectral  regions  and  media,  the  concept  of  the 
modularized  version  was  conceived.  This  version^  was 
designed  to  represent  exactly  the  same  calculations  as  the 
original,  except  for  the  simplification  of  the  program 
structure  into  modules  or  subroutines.  However,  upon  the 
termination  of  that  task  the  authors  proceeded  to  add 
models  for  the  trace  gases,  as  developed  during  the  pre¬ 
sent  scientific  effort. 

7.2  Structure  of  Modularized  Version 

The  basic  design  used  was  that  of  a  main  program 
which  reads  input  data,  computes  total  transmittance  and 
radiance  and  generates  outputs,  and  a  series  of  subroutines 
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which  select  individual  models  and  compute  individual  trans- 
mittances  and  absorber  amounts.  This  is  shown  in  Fig.  4. 
The  main  operational  flow  chart  follows  in  Fig.  5.  Ex¬ 
cluding  the  four  subroutines  for  the  trace  gases,  the 
modularized  version  breaks  down  the  original  into  one 
program  with  11  subroutines.  The  flow  chart  for  subroutine 
ABSORB  is  shown  in  Fig.  6.  This  subroutine  computes  the 
equivalent  absorber  amount  for  all  of  the  attenuation 
models  according  to  Eq.(4),  which  in  terms  of  the  meteo¬ 
rological  variables  becomes 

p (7  \  n  T  m 

w  -  J  (  pf)  <Tify>  du  (56) 

Figure  7  gives  details  of  the  Transmittance/Radiance  Loop 
of  program  Main.  It  is  worth  noting  that  the  modularized 
version  of  Lowtran  being  done  by  AFGL  separates  this  loop 
into  a  subprogram.  The  modularization  discussed  in  this 
text  leaves  the  loop  as  part  of  the  main  program,  but 
extracts  individual  subroutines  for  the  calculation  of  the 
equivalent  absorber  amount,  the  frequency  selection,  and 
the  attenuation  models. 

The  flow  chart  for  FREQSL  subroutine  is  shown  in 
Fig.  8.  This  subroutine  is  designed  to  simplify  the  pro¬ 
cess  of  arriving  at  the  individual  models  effective  at 
the  frequency  of  interest.  It  should  also  assist  the 
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Conceptual  flow  chart  of  modularized  Lowtran 
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Subroutine  ABSORB 


©  © 


Fig,  6 ,  Flow  chart  for  subroutine 

ABSORB  t  computing  equivalent 
path  quantities 


X-HMIN 
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© 


© 


Fig  .  6 


(cont'd) 


Compute  partial  radiance  sums  from  layers 
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Fig.  8. 


Flow  chart  for  subroutine  FREQSL  and  thirteen  transmittance 

model  subroutines. 
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Go  to  40 


Fig.  8 
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Fig,  8 . 
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If  (I.GE.278  and  I.LE,282)^ 
or  (I.GE.326  and  I,LE.346) 
or  (I,GE,480  and  I.LE,498) 
or  (I.GE.511  and  I,LE.S3l)^ 

- > 

It  (I.GE.283  and  I.LE.325) 
_ / 

CALL 

LUAP 

DIVAD 
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NHOJ 

EKUL 

Go 


CALI 
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NHOJ 

EKUL 

ARZE 


If  (I.GE.347  and  I.LE.420) > 
or  (I.GE.439  and  I.LE.479)  , 

If  (I.GE.421  and  I.I.E.438)Ss 
_ ^ 

J  If  (I.GE.531  and  I.LE.585) 

y 
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EVETS 
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NHOJ 
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NHOJ 

EKUL 

SUSEJ 
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SEMAJ 
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"TZ 


?o 


Go 


*> 

to 


Subroutine  SUSEJ 
SO, 

P*C12+LogW (12) 
Tx(12)-EXP[-10**(51  + 
52*p)  ] 


40 


Subroutine  ARZE 
NO 

P-C13+LogW(13) 
Tx(13)-EXP[-lO**(FNl  + 
FN2*p) ] 


40 


Subroutine  SUTIT 
NH3  ' 

P*C.l4+LogW(14) 
Tx(14)«=EXP[“10**(EHl  + 
FH2*p) 1 


Go  to  40 


Subroutine  BOJ 

no2 

P=C15+LogW(15) 
Tx(15)«EXP[-10**(01  + 

02*p) ] 


Go  to  40 


Go  to  40 


Fig,  8 . 
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user  who  desires  to  replace  or  add  models  to  the  program 
and  It  should  reduce  the  overall  computational  time.  This 
subroutine  is  based  on  Tables  2  and  8, 

7.3  Models  for  H»0  Vapor,  Infrared  0_  and  the  Uniformly- 
Mixed  Gases 

As  indicated  abov§,  all  the  attenuation  models  were 

extracted  from  the  main  program  and  placed  into  subroutines. 

The  models  were  left  basically  in  the  same  structural  form 

except  for  the  models  for  HNO^  and  1^0  vapor,  infrared  0^ 

and  the  uniformly-mixed  gases.  The  change  in  the  first 

was  to  arrange  it  along  the  same  form  as  in  the  other 

models  originally  available  in  Lowtran,  That  is,  the 

spectral  parameters  were  extracted  from  the  subroutine 

and  read  at  the  beginning  of  program  MAIN.  The  changes 

in  the  latter  three  gases  (i.e.  ^0  vapor,  0^  in  the 

infrared  and  the  uniformly-mixed  gases)  were  based  on  a 

2 

previous  work  by  Pierluissi  et  al.  on  the  representation 
of  the  tabulated  transmission  functions  by  analytical 
functions.  The  other  principal  change  consisted  of  adding 
models  for  the  trace  gases  SO^*  NO,  NO2  and  NH^. 

To  arrive  at  the  analytical  function  for  modeling 
1^0  vapor  and  the  uniformly-mixed  gases  the  double  expo¬ 
nential  expression 

T  -  exp(-10ao+alX) 

where  x  is  as  in  Eq .  9  and  aQ  and  a^  are  absorber 


(57) 
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constants,  was  curve-fitted  to  the  134  values  of  t  and  x 

tabulated  in  Lowtran.  The  values  found  are  a  =  -1.14619 

o 

and  a^  =  0.55013,  and  it  reproduced  the  tabulated  trans¬ 
mittance  with  a  standard  deviation  of  0.005.  For  0^  the 
function  adopted  is  given  by 


where  a  =  -3.08019,  a,  *  2.11127,  and  the  tabulated  data 
is  reproduced  with  a  standard  deviation  of  0.007.  Note 
in  each  one  of  these  functions  that  the  134  tabulated 
values  are  replaced  with  two  and,  hence,  their  adoption 
reduces  the  computer  storage  requirements.  Also,  they 
inherently  offer  exponential  interpolation  while  with 
the  present  tabulation  linear  interpolation  is  being 
used.  Finally,  there  is  no  need  for  the  small  optical 
thickness  (i.e.  0.999  <  T  <  1)  correction  inserted  in 

Lowtran  4,  as  required  by  its  radiance  calculational 
scheme . 


7.4  Models  for  Trace  Gases  S02>  NO,  NO2,  and  NH^ 

Absorption  ^y  the  trace  gases  was  incorporated  in 
Lowtran  using  a  somewhat  similar  procedure.  Empirical 
transmission  functions  were  first  obtained  from  a  com¬ 
puterized  procedure  which  replaced  the  classical  manual 
graphical  techniques.  The  procedure  is  explained  in 
Chapter  III  of  this  report  and  has  been  proposed  to  the 
scientific  community^. 
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Instead  of  either  representing  the  transmission  function 
by  a  table  or  by  a  single  function,  it  was  divided  into 
nine  segments  for  each  absorber.  The  individual  curve 
segments  are  summarized  in  Table  5,  each  one  being  repre¬ 
sented  by  the  function 


T 


a  +  ax 

exp  (r-10  ) 


(59) 


For  each  absorber  x  is  computed  with  Eqs.  (8)  through  (11) 
and  the  relation 

U  =  0.772  x  10  ^  ppm  p  Z  (60) 

d 


where  ppm  is  the  parts  per  million  by  volume,  is  the 
air  density  in  gm/m  and  Z  is  the  range  in  kilometers. 

Table  6  lists  the  ppm  and  temperature  and  pressure  exponents 
used  in  the  modularized  code  for  the  individual  trace 
gases.  The  ppm  values  are  read  as  input  through  a  separate 
card  which  may  be  easily  changed  according  to  the  needs  of 
the  user.  The  constants  C*  are  tabulated  in  Table  7.  The 
spectral  coverage  for  each  gas  is  depicted  in  Table  8. 

The  models  are  for  a  resolution  of  20  cm  ^  and  are  defined 
at  5  cm  ^  through  their  spectral  regions  of  effectiveness. 
Their  mean  standard  deviation  in  fitting  the  original  line- 
by-line  data  is  about  0.008.  Figure  9  depicts  the  trans¬ 
mission  functions  for  the  four  trace  gases  considered. 


****•*'*>*4 


gments  in  the  empirical  transmission  function 
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TRACE 

GAS 

SPECTRAL 
REGION 
(cm--*- ) 

PRESSURE 

EXPONENT 

n 

TEMPERATURE 

EXPONENT 

m 

PARTS  PER 
MILLION  BY 
VOLUME  ppm 

so2 

440-  615 
1,055-1,250 
1,310-1,410 

0.07122 

0.06159 

0.221 

NO 

1,760-1,970 

0.90098 

1.01192 

0.250 

no2 

1 

655-  880 
1,540-1,670 
2,840-2,895 

0.18066 

0.20911 

0.090 

nh3 

670-1,230 

0.52125 

-0.60438 

0.200 

Table  6 

.  Absorber  parameters 

with  the  models  for 

in  Modularized  Lowtran  usi 
the  trace  gases. 

WAVENUMBER 
(cm-1 ) 

n 

WAVENUMBER 

(cm-1) 

B 

WAVENUMBER 

(cm-1) 

B 

440 

-2.987 

1070 

-1.653 

1320 

-1.237 

445 

-2.330 

1075 

-1.443 

1325 

-0.494 

450 

-1.791 

1080 

-1.252 

1330 

0.139 

455 

-1.370 

1085 

-1.080 

1335 

0.613 

460 

-1.041 

1090 

-0.926 

1340 

0.899 

465 

-0.795 

1095 

-0.787 

1345 

1.043 

470 

-0.613 

1100 

-0.661 

1350 

1.090 

475 

-0.469 

1105 

-0.544 

1355 

1.097 

480 

-0.346 

1110 

-0.434 

1360 

1.104 

485 

-0.233 

1115 

-0.329 

1365 

1.093 

490 

-0.126 

1120 

-0.230 

1370 

1.118 

495 

-0.037 

1125 

-0.139 

1375 

1.088 

500 

o 

o 

1130 

-0.073 

1380 

0.926 

505 

-0.008 

1135 

-0.047 

1385 

0.534 

510 

-0.052 

1140 

-0.057 

1390 

-0.067 

515 

-0.102 

1145 

-0.083 

1395 

-0.804 

520 

-0.102 

1150 

-0.098 

1400 

-0.768 

525 

-0.044 

1155 

-0.071 

1405 

-1.687 

530 

0.013 

1160 

-0.020 

1410 

-2.469 

535 

0.039 

1165 

0.014 

2450 

-3.669 

540 

0.014 

1170 

0.011 

2455 

-2.855 

545 

-0.056 

1175 

-0.040 

2460 

-2 . 131 

550 

-0.141 

1180 

-0.123 

2465 

-1.528 

The  spectral  coefficient  C'(v)  for  SOj* 


Table  7a 
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WAVENUMBER 

(cm*^) 

C' 

WAVENUMBER 

(cm"1) 

C’ 

WAVENUMBER 

(cm-1) 

C' 

555 

-0.221 

1185 

-0.213 

2470 

-1.076 

560 

-0.294 

1190 

-0.301 

2475 

-0.805 

565 

-0.366 

1195 

-0.388 

2480 

-0.647 

570 

-0.442 

1200 

-0.481 

2485 

-0.571 

575 

-0.529 

1205 

-0.586 

2490 

-0.549 

580 

-0.635 

1210 

-0.707 

2495 

-0.539 

585 

-0.766 

1215 

-0.843 

2500 

-0.536 

590 

-0.934 

1220 

-0.996 

2505 

-0.517 

595 

-1.157 

1225 

-1.165 

2510 

-0.528 

600 

-1.457 

1230 

-1.351 

2515 

-0.691 

605 

-1.862 

1235 

-1.554 

2520 

-1.073 

610 

-2.420 

1240 

-1.777 

2525 

-1.673 

615 

-3.094 

1245 

-2.033 

2530 

-2.414 

1055 

-2.604 

1250 

-2.369 

2535 

-2 . 207 

1060 

-2.156 

1310 

-3.010 

1065 

/ 

-1.884 

1315 

-2.080 

Table  7a,  (Continued) 
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efficient  C'(v)  for  NO 
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WAVENUMBER 

C ' 

WAVENUMBER 

C ' 

WAVENUMBER 

C' 

655 

-0.844 

800 

m 

1,600 

2.616 

660 

-0.760 

805 

-0.286 

1,605 

2.616 

665 

-0.676 

810 

-0.315 

1,610 

2.606 

670 

-0.608 

815 

-0.334 

1,615 

2.608 

'675 

-0.543 

820 

-0.352 

1,620 

2.643 

680 

-0.496 

825 

-0.366 

1,625 

2.682 

685 

-0.450 

830 

-0.396 

1,630 

2.672 

690 

-0.414 

835 

-0.423 

1,635 

2.576 

695 

-0.383 

840 

-0.459 

1,640 

2.350 

700 

-0.326 

845 

-0.498 

1,645 

1.955 

705 

-0.289 

850 

-0.541 

1,650 

1.346 

710 

-0.217 

855 

-0.586 

1,655 

0.596 

715 

-0.140 

860 

-0.630 

1,660 

-0.258 

720 

-0.097 

865 

-0.676 

1,665 

-1.214 

725 

-0.034 

870 

-0.720 

1,670 

-1.951 

730 

-0.031 

875 

-0.766 

2,840 

-1.220 

735 

-0.082 

880 

-0.809 

2,845 

-0.644 

740 

-0.139 

1,540 

-2.428 

2,850 

-0.253 

745 

-0.216 

1,545 

-1.494 

2,855 

0.052 

750 

-0.249 

1,550 

-0.647 

2,860 

0.326 

755 

-0.207 

1,555 

0.122 

2,865 

0.574 

760 

-0.117 

1,560 

0.756 

2,870 

0.792 

The  spectral  coefficient  C?(v)  for  N02> 


Table  7c 
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■H 

WAVENUMBER 

C' 

WAVENUMBER 

C' 

765 

-0.047 

1,565 

1.230 

2,875 

0.978 

770 

0.000 

1,570 

1.568 

2,880 

1.122 

775 

0.009 

1,575 

1.855 

2,885 

1.216 

780 

-0.046 

1,580 

2.104 

2,890 

1.252 

785 

-0.100 

1,585 

2.310 

2,895 

1.249 

790 

-0.148 

1,590 

2.469 

795 

-0.214 

1,595 

2.573 

Table  7c 


(Continued) 


WAVENUMBER 

C  ' 

WAVENUMBER 

mm 

WAVENUMBER 

C' 

690 

-2.603 

875 

-1.124 

1,060 

-0.589 

695 

-2.456 

880 

-1.155 

1,065 

-0.565 

700 

-2.290 

885 

-1.161 

1,070 

-0.537 

705 

-2.128 

890 

-1.143 

1,075 

-0.510 

710 

-1.980 

895 

-1.139 

1,080 

-0.512 

715 

-2.225 

900 

-1.117 

1,085 

-0.528 

720 

-1.823 

905 

-1.107 

1,090 

-0.575 

725 

-1.744 

910 

-0.844 

1,095 

-0.625 

730 

-1.674 

915 

-0.558 

1,100 

-0.668 

735 

-1.577 

920 

-0.238 

1,105 

-0.694 

740 

-1.481 

925 

-0.042 

1,110 

-0.717 

745 

-1.372 

930 

-0. 002 

1,115 

-0.740 

750 

-1.284 

935 

-0.157 

1,120 

-0.774 

755 

-1.207 

940 

-0.436 

1,125 

-0.834 

760 

-1.128 

945 

-0.610 

1,130 

-0.905 

765 

-1.061 

950 

-0.548 

1,135 

-0.977 

770 

-1.004 

955 

-0.352 

1,140 

-1.042 

775 

-0.947 

960 

-0.139 

1,145 

-1.133 

780 

-0.886 

965 

-0.095 

1,150 

-1.219 

785 

-0.876 

970 

-0.365 

1,155 

-1.301 

790 

-0.872 

975 

-0.729 

1,160 

-1.383 

795 

-0.869 

980 

-1.048 

1,165 

-1.488 

800 

-0.872 

985 

-1.275 

1,170 

-1.594 

Table  7d.  The  spectral  coefficient  C*(v)  for  NH 
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WAVENUMBER 

Hi 

WAVENUMBER 

mm 

WAVENUMBER 

c* 

805 

-0.848 

990 

-1.257 

1,175 

-1.696 

810 

-0.811 

995 

-1.142 

1,180 

-1.796 

815 

-0.772 

1,000 

-1.053 

1,185 

-1.873 

820 

-0.773 

1,005 

-0.963 

1,190 

-1.936 

825 

-0.793 

1,010 

-0.920 

1,195 

-1.991 

830 

-0.825 

1,015 

-0.944 

1,200 

-2.080 

835 

-0.869 

1,020 

-0.889 

1,205 

-2.183 

840 

-0.894 

1,025 

-0.829 

1,210 

-2.292 

845 

-0.890 

1,030 

-0.736 

1,215 

-2.404 

850 

-0.873 

1,035 

-0.644 

1,220 

-2.529 

855 

-0.868 

1,040 

-0.596 

1,225 

-2.639 

860 

-0.907 

1,045 

-0.569 

1,230 

-2.732 

865 

-0.965 

1,050 

-0.572 

870 

-1.045 

1,055 

-0.590 

_ 

Table  7d 


(Continued) 
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Fig. 


Empirical  transmission  function  for  NO. 
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VIII.  Calculations  and  Results 
8.1  Introduction 

The  procedure  for  the  use  of  the  Modularized  Lowtran 
in  calculations  is  identical  to  that  of  the  original  and, 
hence,  deserves  no  further  explanation.  There  are  some 
input  and  output  alterations  that  deserve  some  explanatory 
remarks.  Changes  in  the  input  format  include: 

1.  Reading  of  the  spectral  constants  for  all  band 
models  at  the  beginning  of  the  main  program 
rather  than  in  the  subroutines. 

2.  Elimination  of  the  transmittance  tables  for 
H^O  vapor,  infrared  0^  and  the  uniformly-mixed 

gases . 

3.  Reading  of  the  spectral  constants  for  the  newly 
added  band  models  for  the  trace  gases. 

4.  Reading  of  the  air  density  profile  for  the  U.S. 
Standard  atmosphere,  and  of  the  ppm  for  the 
calculation  of  the  equivalent  amounts  of  the 
trace  gases. 

5.  Changes  in  the  dimension  statements  to  include 
the  additional  subscripted  variables. 

Changes  in  the  output  format  include: 

1.  Modification  of  the  print  out  of  the  input  data. 

2.  Modification  of  the  output  table  of  computations 
to  include  the  transmittance  for  the  trace  gases. 

It  should  be  stressed,  however,  that  the  code  is  operated 
using  exactly  the  same  four  control  cards  as  in  the  original 


code . 
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8.2  Testing  of  Modularized  Version 

The  first  step  in  the  testing  of  the  modularization 
consisted  of  running  identical  calculations  using  the 
original  code  and  the  modularized  code  before  the  replace¬ 
ment  of  the  transmittance  tables  and  before  the  addition 
of  the  trace  gases.  Numerous  cases  were  considered  during 
this  effort.  A  particular  case  in  which  the  spectral  range 
varied  from  2350  to  2450  cm  *  for  a  path  at  65°  from  a 
height  of  2.5  km  to  a  height  of  8.5  km  and  a  23  km  visual 
range,  is  shown  in  the  Appendix.  This  output  is  identical 
to  the  output  obtained  from  the  original  Lowtran. 

The  second  step  in  the  testing  procedure  consisted 
of  running  calculations  using  the  original  code  and  the 
modularized  version  with  the  transmittance  tables  replaced 
with  the  continuous  functions,  but  before  the  addition 
of  the  trace  gases.  For  this  purpose,  10  frequencies 
were  selected  such  that  different  combinations  of  models 
would  be  effective  in  the  calculation  of  the  total  trans¬ 
mittance.  The  calculations  were  for  a  5  km  path  at  sea 
level  in  a  sub-arctic  winter  atmosphere  with  a  23  km 
visibility.  The  results  are  summarized  in  Table  9.  The 
columns  listed  under  Transmittance  Deviations  represent 
the  differences  between  the  calculations  using  the  tab¬ 
ulated  and  the  continuous  functions.  Note  that  the  average 
total  transmittance  deviation  is  0.0034,  which  is  below 
the  standard  deviation  obtained  in  the  curve  fitting  of 


TRANSMITTANCE  DEVIATIONS 


WAVENUMBER 

(cm-^) 

h2o 

VAPOR 

1 

1 

INFRARED 

°3 

UNIFORMLY- 
MIXED  GASES 

TOTAL 

TRANSMITTANCE 

455 

0.0022 

0.0000 

0.0000 

0.0021 

555 

0.0035 

0.0000 

0.0026 

0.0018 

655 

0.0041 

0.0003 

0.0000 

0.0000 

755 

0.0007 

0.0003 

0.0047 

0.0038 

955 

0.0057 

0.0002 

0.0050 

0.0096 

1155 

0.0026 

0.0003 

0.0050 

0.0058 

1355 

0.0044 

0.0000 

0.0013 

0.0006 

1855 

0.0007 

0.0001 

0.0034 

0.0007 

2455 

0.0015 

0.0000 

0.0054 

0.0045 

3155 

0.0037 

0.0001 

0.0027 

0.0053 

Transmittance  difference  between  calculations 
using  the  tabulation  of  the  transmittance  func¬ 
tions  and  calculations  using  the  continuous 
function  representation  for  a  5  km  path  at  sea 
level  in  a  sub-arctic  winter  atmosphere. 


Table  9 


the  functions  to  the  individual  transmittance  tables.  This 
deviation  amounts  to  an  error  of  about  0.7%  in  the  middle 
of  the  curve-of -growth ,  which  far  exceeds  the  accuracy  of 
Lowtran  (between  10  to  20%).  The  following  are  attractive 
features  of  the  continuous  functions: 

1.  They  inherently  provide  for  continuous  expo¬ 
nential  interpolation  in  transmittance,  which 
is  superior  to  the  linear  interpolation  used 
in  connection  with  the  transmittance  tables. 

2.  They  provide  for  analytical  operations  such  as 
differentiation  and  interpolation  often  needed 
in  radioactive  transfer  problems. 

3.  They  can  be  used  easily  for  curve  fitting  to 
new  transmittance  data  using  computerised 
procedures . 

4.  Their  use  reduces  significantly  the  computer 
storage  requirements  for  the  individual  models. 

5.  They  continuously  provide  for  transmittance 
calculations  for  small  argument  values  where 
0.9999  <  x  <  1,  for  which  range  Lowtran  4 
includes  an  additional  exponential  function. 

It  should  be  pointed  out  that  the  deviations  listed  in 
Table  9,  although  insignificant,  do  not  represent  errors 
solely  attributed  to  the  analytical  functions.  Since  they 
are  smaller  than  the  uncertainties  in  the  original  data 
used  to  develop  the  tabulated  transmittances ,  they  pri¬ 
marily  represent  differences  in  the  calculational  procedures. 
In  fact,  in  the  region  between  the  tabulations  the  use  of 
the  analytical  functions  are  likely  to  provide  more  ac¬ 
curate  results  than  the  use  of  the  original  method  in  Lowtran 
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The  last  effort  in  the  testing  of  the  modularized 
code  consisted  of  calculations  involving  the  newly  added 
trace  gases.  For  this  purpose  ten  frequencies  were  run 
at  which  the  trace  gas  models  are  effective.  The  same 
frequencies  were  run  with  the  modularized  Lowtran  without 
these  models.  The  results  are  summarized  in  Table  10. 

The  table  is  primarily  intended  to  show  the  absorptive 
effects  of  the  trace  gases. 


total  transmittance  include  all  the  attenuators  and  the  tra 
ept  for  the  rightmost  column  which  excludes  the  trace  gases 


11:  (a)  Atmospheric  regions  Included  In  the  data 

calculations 


Model 

P  Onbar) 

T  (#K) 

Standard 

1013 

288,1 

898.6 

281.6 

795.0 

275.1 

701,2 

268.7 

616.6 

262.2 

Tropical 

805.0 

288.0 

Subarctic 

Winter 

1013 

257.1 

(b)  Transmittance  cuts  chosen  from  the  curve  of 
growth 


Coefficients  of 

Absorber  Spectral  Parameter  Cr  (cm-1)  Analytical  Standard 

Parameters  Function  Deviation 


Coefficients  of 

Absorber  Spectral  Parameter  C'  (crn'^)  Analytical  Standard 

Parameters  Function  Deviation 


-0.0866  0.2064-0.0759  -0.3133  -0.5107  -0.7296  -0.8199  -0.7357  0.1179 


Coefficients  of 

Absorber  Spectral  Parameter  C*  (cm~  )  Analytical  Standard 

Parameters  Function  Deviation 


Table  12c,  Band  model  parameters  for  NO, 


-0.0722  0.6120  0.4010  0.0559  -0,2290  -0.4530-0.5734  -0.4794  0.0785 


4 

* 
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8.3  Band  Model  Development 

Two  sets  of  curves  of  growth  data  for  each  major 
absorption  band  for  four  trace  gases  SOj  NO,  NOj*  and  NH^ 
were  generated  by  the  line-by-line  calculation  from  the 
AFGL  trace  gas  parameter  tape.  One  of  them  consists  of 
12-cut  data  for  several  layers  of  atmosphere  and  the  other 
consists  of  65-cvit  data  for  the  standard  atmosphere  only. 
Considering  the  wide  range  of  applications,  we  included 
not  only  the  standard  atmospheric  conditions  but  also  one 
condition  each  from  the  tropical  and  subarctic  winter 
climates.  They  are  listed  in  Table  11  together  with  the 
12  chosen  transmittance  values.  The  major  absorption 
bands  for  the  four  trace  gases  are  given  in  Table  12 
together  with  the  corresponding  computed  C'  values. 

Ten  middle  cuts  were  chosen  from  the  12-cut  data 
and  used  in  both  ADSET  and  SIMMIN  for  the  computation  of 
the  band  model  parameters  and  the  standard  transmission 
function.  Depending  on  the  number  of  major  absorption 
bands,  the  total  numbers  of  data  used  differ  but  are  in 
the  range  of  60-210.  The  65-cut  data  was  used  in  ADSET 
for  the  piecewise  interpolation  to  compute  piecewise 
analytical  transmission  functions. 

The  ADSET  computations  were  done  first.  The  ob¬ 
tained  band  model  parameter  values  n,  m,  and  and  nine 
sets  of  coefficients  a^,  a^,  and  a^  are  tabulated  in 
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Table  12  for  NO,  NO2,  and  NHj  in  this  order.  The  cor¬ 

responding  standard  deviations  are  alsb  listed  in  these 
tables . 

We  also  have  generated  standard  atmospheric  condition 
data  for  non-major  bands  of  each  trace  gases.  These  data 
were  used  to  evaluate  non-major  C'(v)  values.  The  computed 
C.'(v)  values  were  listed  in  Table  5.  As  we  have  discussed, 
these  C'(v)  values  and  the  band  model  parameters  together 
with  the  first  order  piecewlse-analy tical  standard  trans¬ 
mission  function  were  implemented  in  the  modularized 
Lowtran. 

We  recall  that  the  S1MMIN  computation  is  a  recursive 
one  and  we  need  a  set  of  initial  guesses  of  the  parameter 
values  to  start  the  computation.  For  the  band  model 
parameters  n,  m,  and  C|,  we  used  the  values  computed  by 
ADSET.  For  a^  and  a^,  the  respective  averages  of  the  first 
order  piecewise  interpolation  results  of  ADSET  were  used. 
Finally,  a^  was  set  to  be  zero.  We  note  that  our  initial 
guesses  are  fairly  accurate,  since  these  values  were 
optimal  or  optimal  in  average  for  ADSET  computation.  A 
small  number  £  which  was  used  for  the  check  of  convergence 
was  chosen  to  be  10  Since  the  parameter  values  are 

expected  to  be  in  the  range  -10-10,  £  -  10  ®  gives  the 
limit  of  numerical  accuracy  of  numbers  in  the  computer. 

The  SIMMIN  results  are  also  listed  in  Table  12. 


rat 


Typical  curve-fits  by  piecewise  analytical  standard 


transmission  functions. to  actual  data  are  shown  in  Fig. 

10  for  SO2  at  500  wavenumber.  The  corresponding  analytical 
standard  transmission  function  are  also  compared  to  the 
data  in  Fig.  10.  In  all  of  the  three  graphs  in  this 
Figure,  the  65-cut  data  were  also  plotted  to  show  the 
fitness  of  the  standard  curves. 


The 

computation  was 

repeated 

using 

two 

smaller 

data 

sets 

with  6 

and  4  cuts  only. 

The  chosen 

cuts  were 

(0.95 

,  0.9 

,  0.8, 

0.6,  0.4, 

and  0.1) 

for  6 

cut 

data  and 

(0.95 

,  0.9 

,  0.6, 

and  0,2) 

for  4  cut 

data . 

The 

derived 

band  model  parameter  values  were  similar  to  those  in 
Table  12  and,  hence,  were  not  repeated  here.  Instead, 
the  corresponding  standard  deviations  were  listed  and 
compared  with  the  10  cut  cases  in  Table  13. 


STANDARD  DEVIATIONS  IN 


Table  13.  Comparison  of  standard  deviations  in  t.  The  two  rows,  ADSET  1  and  2  are,  re 
spectively,  for  the  piecewise  analytical  transmission  functions  and 
linear  and  quadratic  exponents. 
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IX  Discussion  and  Conclusions 

9.1  Introduction 

The  modularized  version  presented  here  is  funda¬ 
mentally  the  same  Lowtran  code  except  for  the  separation 
of  its  computation  structure  into  separate  modules  or 
subroutines.  Although  if  is  based  on  the  4th  version,  it 
can  be  adapted  with  little  modification  to  any  future 
versions,  such  as  the  Sth  version  now  in  progress.  In 
fact,  this  latter  version  already  has  been  structured  by 
AFGL  such  that  the  emission/radiance  loop  is  in  a  subroutine. 
The  modularized  code  presented  here  breaks  down  that  loop 
into  a  frequency  selection  subroutine,  an  equivalent 
absorber  amount  subroutine  and  separate  subroutines  for 
each  one  of  the  attenuation  codes.  The  use  of  modules  in 
a  complex  code  such  as  Lowtran  has  numerous  advantages, 
among  which  the  amenability  for  updating  by  individual 
users  to  suit  their  specific  needs  is  at  the  top  of  the 
list.  In  the  ever  changing  field  of  modeling  it  is  highly 
desirable  to  be  able  to  easily  modify  the  code  for  changes 
in  the  spectral  coverage,  the  spectral  resolution,  the 
absorber  concentrations  in  abnormal  environments,  the 
original  transmission  data  used  in  the  development  and  in 
the  models  used  for  the  individual  attenuators.  The  mod¬ 
ularized  version  presented  here,  although  is  not  the  final 
answer  to  all  conceivable  needs,  it  is  a  first  basic  step 
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in  that  direction.  Practicing  this  predicament,  the 
authors  added  transmission  models  for  the  trace  gases  to 
the  code. 

9.2  Changes 

The  following  are  the  basic  changes  introduced  in 
The  Modularized  Lowtran: 

1.  The  original  main  program  was  separated  into  a 
central  program  and  subroutines  for  the  absorber 
amount  and  the  Individual  attenuation  models. 

2.  In  the  interest  of  efficiency  and  clarity,  a 
new  subroutine  FGQSL  was  added  for  the  selection 
of  the  attenuation  model  effective  at  the  given 
frequency . 

3.  The  subroutine  HNO^  was  re-structured  to  the 
form  of  the  other  previously  incorporated  sub¬ 
routines  in  Lowtran. 

4.  Continuous  analytical  models  were  provided  to 

replace  the  transmittance  curves  for  H^O  vapor, 
Og  and  the  uniformly-mixed  gases.  * 

5.  New  subroutines  for  the  trace  gases  SC^,  NO,  NOj 
and  NH^  were  added. 

A  copy  of  the  modularized  version  is  found  in  the  Appendix. 


9.3  Model  Development 

The  values  of  band  model  parameters  n  and  m  and 
spectral  parameters  obtained  by  ADSET  and  SIMMTK  agreed 
very  well.  Turthermore,  as  It  was  shown  In  Table  13*  the 
standard  deviations  corresponding  to  different  cases 
followed  a  same  pattern  for  the  ADSET  and  SIMMXN  results. 
This  consistency  proves  the  validity  of  both  methods. 

In  general,  the  SIMMIN  and  ADSET  computations  re¬ 
sulted  in  similar  standard  deviations.  It  was  expected 
that  the  ADSET  computation  should  result  in  lower  standard 
deviations  since  it  contained  more  parameters  to  adjust. 
However,  for  a  half  of  the  cases,  the  SIMMIN  code  produced 
lower  standard  deviations.  This  is  due  to  the  large 
computational  error  for  the  ADSET  computations  in  solving 
the  normal  equation  AX  ■  B.  When  the  condition  number  of 
the  coefficient  matrix  A  becomes  large  (i.e.,  A  becomes 
close  to  be  singular),  the  computational  error  becomes  so 
large  that  it  can  exceed  the  directly  minimized  error  of 
the  SIMMIN  computation. 

We  note  that  this  reversal  occurred  for  all  four 
cut  data  cases.  This  suggests  that  the  advantage  for  ADSET 
of  having  more  parameters  to  be  adjusted  is  not  significant 
for  these  cases.  Hence,  we  recommend  the  use  of  SIMMIN 
if  the  available  data  contains  less  than  five  or  six  cuts. 

A  comparison  of  the  standard  deviations  for  two 


piecewise  interpolation  results  in  the  ADSET  computation 
showed  no  significant  difference.  Furthermore,  the  results 
with  the  second  method  using  quadratic  form  of  x  on  the 
exponent  of  the  double  exponential  function  were  'bumpy* 
for  some  cases.  Since  the  nature  of  the  transmittance 
does  not  predict  this  behavior,  we  conclude  that  the  first 
method  using  linear  function  of  x  is  accurate  enough  to  be 
used  in  the  actual  application. 

The  standard  deviations  were  much  higher  for 
cases  than  the  cases  for  the  rest  of  absorbers.  By  in¬ 
specting  each  curve  of  growth  in  detail,  it  was  found  that 
this  was  mainly  due  to  the  difference  in  the  steepness  of 
the  curves  of  growth  for  three  absorption  bands.  This 
difference  cannot  be  compensated  by  values  since  they 
only  shift  the  curves  of  growth  linearly.  In  fact,  within 
the  current  band  model  structure,  it  is  Impossible  to 
compensate  this  difference.  Hence,  it  may  be  necessary  to 
modify  some  of  the  basic  assumptions  regarding  the  band 
model  structure,  if  lower  standard  deviations  are  required. 

As  a  side-effect  of  this  discrepancy  in  the  tangent 
of  curves  of  growth,  the  SIMMIN  computation  took  far  more 
time  for  NO^  cases  than  the  rest.  Most  of  the  computations 
of  ADSET  were  completed  by  26-36  CPU  seconds.  The  fluc¬ 
tuations  in  the  computation  time  were  very  small.  On  the 
other  hand,  the  SIMMIN  computation  time  varied  from  14 
seconds  to  270  seconds.  NO^  cases  consumed  about  200-270 


seconds;  which  were  about  four  times  as  much  as  that  for 
the  other  cases.  This  is  because  the  minimizing  point  in 
the  parameter  space  is  not  well  defined  for  N<>2  cases.  In 
other  words,  the  error  surface  in  thi  parameter  space  has 
a  very  shallow  bottom  so  that  the  updating  step  cannot 
produce  large  enough  changes  in  the  parameter  guesses  in 
order  to  have  a  rapid  convergence. 

Thus,  it  was  found  that  the  accuracy  of  the  computed 
results  and  the  time  of  execution  depend  heavily  on  the 
actual  data.  Hence,  it  is  very  important  to  give  enough 
consideration  for  the  data  structure.  This  will  be  dis- 


98 


9.4  Data  Structure 

As  it  was  expressed  earlier,  we  assumed  that  the 
number  of  layers  (=  the  number  of  data  points)  in  each  cut 
is  the  same  for  an  absorption  band.  This  was  done  for  the 
sake  of  easier  coding  in  data  handling.  However,  this 
assumption  need  not  be  valid.  Especially  in  weaker  absorp¬ 
tion  bands,  it  is  required  to  use  very  large  range  values 
to  have  high  enough  equivalent  absorber  amounts  in  order 
to  realize  lower  transmittances .  In  some  cases  the  range 
becomes  enormous  (in  the  order  of  the  radius  of  the  earth) 
so  that  the  corresponding  data  no  longer  possess  physical 
significance.  The  ADSET  code  has  a  criterion  that  if  the 
logarithm  of  the  equivalent  absorber,  log  W,  exceeds  a 
certain  critical  value,  then  the  corresponding  data  will 
be  set  aside  and  will  not  be  used  in  the  later  computation. 
The  critical  value  was  set  to  be  2  for  the  actual  computa¬ 
tion,  which  corresponds  approximately  to  a  vertical  path 
through  the  atmosphere. 

In  connection  with  this,  if  data  are  not  available 
at  some  layers,  then  the  data  values  are  set  at  0  to  flag 
the  nonavailability  of  data.  ADSET  can  also  detect  this 
and  will  ignore  the  data. 

A  caution  must  be  executed  in  choosing  combinations 
of  pressures  and  temperatures,  i.e.,  atmospheric  conditions 
If  a  data  set  contains  either  the  standard  pressure  or  the 


standard  temperature  or  both  only,  then  both  ADSET  and 
SIMMIN  fall  because  of  the  fact  that  the  coefficient  of  n 
or  m  or  both  In  Eq.(12)  becomes  zero,  since 

P 

log  Cp5-)  -  log  1  -  0,  (61) 

o 

T 

log  GpO  -  log  1  -  0.  (62) 

o 

For  this  case,  the  coefficient  matrix  A  of  the  normal 
equation  In  ADSET  becomes  singular  and  the  gradient  cor¬ 
responding  to  n  or  m  or  both  in  SIMMIN  becomes  zero  all  the 
time.  Hence,  the  normal  equation  cannot  be  solved  In  ADSET 
and  the  initial  guess  of  n  or  m  or  both  cannot  be  changed 
In  SIMMIN. 

Another  consideration  which  should  be  pointed  out 
is  to  Include  different  climate  conditions.  The  standard 
climate  condition  for  several  layers  of  atmosphere  contains 
sequence  of  pressures  and  temperatures  both  of  which  are 
monotone  decreasing.  Therefore,  if  only  these  conditions 
are  used,  then  it  is  very  difficult  to  distinguish  the 
cause  of  changes  in  the  transmittance  due  to  the  changes 
in  pressure  and  in  temperature.  This  leads  to  the  shallow 
bottom  of  the  error  surface  and  hence,  large  computational 
error  results  in  ADSET  caused  by  the  large  condition 
number  of  the  coefficient  matrix  and  slow  convergence  in 
SIMMIN  due  to  the  small  gradient.  In  the  actual  computation, 
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we  included  not  only  the  standard  climate  conditions  but 
also  one  condition  each  from  the  tropical  and  subarctic 
winter  climates  in  consideration  of  wide  applicability  of 
the  results.  Numerically  speaking,  this  also  resulted  in 
making  the  regression  problem  well-posed  by  breaking  the 
monotonousness  of  the  pressure  and  temperature  combinations 
of  the  standard  conditions.  In  fact,  several  computations 
were  done  for  ADSET  and  SIMMIN  with  standard  condition  data 
only.  SIMMIN  took  10-45  minutes  of  CPU  time  to  converge 
if  it  were  convergent  and  ADSET  resulted  in  a  set  of  absurd 
values  for  n  and  m.  Thus,  the  importance  of  the  numerical 
consideration,  which  is  Ignored  in  many  cases,  is  clearly 
indicated.  The  proper  care  should  be  taken  when  selecting 
controllable  data  values. 
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ADSET  V  1-3 

EVALUATION  OF  ABSORBER  PARAMETERS  AND  STANDARD  EMPIRICAL  AND 
PIECEWISE- ANALYTICAL  TRANSMISSION  FUNCTIONS 

THIS  CODE  USES  THE  SUBROUTINE  SIMQ  IN  SSP  LIBRARY 

THIS  CODE  CONSISTS  OF 

MAIN:  COMPUTATIONS  OF  BAND  PARAMETERS  N,M,C' 

AND  EMPIRICAL  STANDARD  TRANSMISSION  FUNCTION 
NMBC :  COMPUTATION  OF  NON-MAJOR  BANDS’  C-VALUES 

INTPL1:  COMPUTATION  OF  THE  PIECEWISE-ANALYT ICAL 

STANDARD  TRANSMISSION  FUNCTION  GIVEN  BY 
TAU  =  E  X P ( —  1 0 *  * ( A 1 +A2*X  ) ) 

INTPL2 :  COMPUTATION  OF  THE  PIECEWISE-ANALYT ICAL 

STANDARD  TRANSMISSION  FUNCTION  GIVEN  BY 
TAU  =  EXP(-10*«(A1+A2«X+A3*X**2) ) 

SDTAU :  COMPUTATION  OF  THE  ERROR  STANDARD  DERIVATIONS 

BETWEEN  PIECEWISE  ANALYTICAL  STANDARD 
TRANSMISSION  FUNCTION  AND  THE  ORIGINAL  DATA 
USED  IN  THE  MAIN  PROGRAMME 

DATA  SFT-IIP 

1.  1ST  CARD  TITLE  IN  20A4 

(ABSORBER  TYPE, ETC) 

2.  2ND  CARD  FOUR  CONTROL  NUMBERS  IN  415 

MAXRPT  MAXIMUN  NUMBER  OF  REPETITION  OF  . 

THE  COMPUTATION  IN  MAIN 
INDX(1)  SUBROUTINE  NMBC  IS  CALLED  IF  >  0 
INDX(2 )  SUBROUTINE  INTPL1  IS  CALLED  IF  >  0 
INDXC3)  SUBROUTINE  INTPL2  IS  CALLED  IF  >  0 

3.  DATA  SET  FOR  MAIN 

CONSISTS  OF  SEVERAL  SUBSETS  (MAX. 6)  OF  DATA 
CORRESPONDING  TO  INDIVIDUAL  MAJOR  BANDS 
FIRST  CARD  FOR  EACH  SUBSET  IS  A  CONTROL  CARD  WHICH 
CONTAINS  BAND#,  WAVE  NUMBER,  #  OF  CUTS  AND  #  OF  LEVELS 
IN  THIS  ORDER  BY  THE  FORMAT  ( 15 , F 1 0 . 3 , 215 ) 

.  REFER  TO  THE  READ(5,105)  AND  FORMAT  105  FOR  THE 
.  CONTENTS  OF  EACH  CARD. 

END  OF  DATA  IS  DEFINED  BY  A  BLANK  CONTROL  CARD 

4.  DATA  SET  FOR  NMBC 

INDIVIDUAL  DATA  FORMAT  -  SAME  AS  FOR  MAIN 
END  OF  EACH  DATA  SET  FOR  A  BAND  IS  MARKED  BY 
A  BLANK  CARD 

END  OF  ALL  DATA  IS  MARKED  BY  -1  (12)  IN  ADDITION  TO 
A  BLANK  CARD 

IF  NO  DATA  BUT  A  BLANK  CARD  IS  SUPPLIED,  THEN 
THIS  SUBROUTINE  IS  SKIPPED. 

5.  DATA  SET  FOR  INTPL1 

FORMATION  OF  THE  DATA  IS  THE  SAME  AS  THAT  FOR  MAIN 
DATA  (IF  SUPPLIED)  WILL  BE  USED  FOR  S.D.  COMPUTATION 
ONLY. 

IF  NO  DATA  BUT  A  BLANK  CARD  IS  SUPPLIED,  THEN  THE 
S.D.  COMPUTATION  IS  SKIPPED. 

6.  DATA  SET  FOR  INTPL2 
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FORMATION  OF  THE  DATA  IS  THE  SAME  AS  THAT  FOR  MAIN 
IF  NO  DATA  BUT  A  BLANK  CARD  IS  SUPPLIED,  THEN 
THIS  SUBROUTINE  IS  SKIPPED. 

NOTE:  DATA  SET  MUST  HAVE  A  CUT  STRUCTURE  SUCH  THAT  EQUAL 

TRANSMITTANCE  DATA  ARE  GROUPED  TOGETHER  AND  THESE  GROUPS 
ARE  QUEUED  IN  THE  DECENDING  ORDER  IN  TAU.  THE  QUEUING  OF 
THE  LEVELS  WITHIN  EVERY  GROUP  MUST  BE  THE  SAME. 

DIMENSION  V(19) ,A(19, 19) ,X(361 ) ,B(19) ,RI(6, 12, 10) 

DIMENSION  P( 1 0 ) , WWW( 12,10), STANDV (12), TSD(6 ) , NDATA(6 ) , INDX( 3 ) 
COMMON  /PARM  1  /  TSTD( 1 2 ) , PW( 1 2 ) , WN (6 ) , CSTD (6 ) , NCUT , NC , NAME (20 ) , 

*  AN,AM,CF,IC0NST(6) ,NEL 

COMMON  /PARM2/  PRES ( 6 , 1 2 , 1 0 ) , TEMP ( 6 , 1 2 , 1 0 ) , UGAS(6 , 1 2 , 1 0 ) , 

*  TAU (6 , 12), NTC ( 6 ) , NLV (6 ) 

CF=1 .0 

L00PCT=1 
WCR IT=2 . 

M=0 

READ(5 , 100 )  (NAME(I) ,1=1 ,20) 

100  FORMAT (20A4) 

READC5, 101 )  MAXRPT, (INDX(I) ,1=1 ,3) 

101  FORMAT (415 ) 

COMPUTATION  OF  ABSORBER  PARAMETERS  N,  M  4  C-VALUES  IS  REPEATED 
MAXRPT  TIMES,  WHERE  1  <  MAXRPT  <  10  IS  READ  IN  BY  15  FORMAT 
(SUGGESTED  VALUE  IS  5) 

DATA  READ-IN  ROUTINE 

1000  CONTINUE 

READ(5 ,102)  IC , W , JM , KM 

102  FORMAT(I5,F10. 3,215) 

IF(IC.LE.O)  GO  TO  2000 
IF(M.GT.O)  GO  TO  10 

CALL  DATE  (MONTH , IDAY , IYEAR) 

WRITE (6, 1 1 1 ) MONTH, IDAY, IYEAR 
111  FORMAT ( 1H1 ,T60, 14, •  /  ’,12,'  /  ’,12,//) 

WRITE (6 , 200 )  (NAME(I) ,1=1 ,20) 

200  FORMAT ( 1H  ,T25,20A4) 

GO  TO  11 

10  CONTINUE 
WRITE(6 ,201  ) 

201  FORMAT ( 1H  1 ) 

11  CONTINUE 
M=M+1 

WN ( M) =W 
NTC ( M) = JM 
NLV ( M) =KM 

WRITE (6 , 202 )  M,WN(M) ,NTC(M) ,NLV(M) 

202  FORMAT ( 1  HO , T 15 ,’ ***  BAND’, 13,'  (WAVE  NUMBER  =’ ,F10 . 3 ,’ )  **,M, 

*  III ,T20 , ’TOTAL  »  OF  CUTS  s’ ,13,//, T20, ’TOTAL  #  OF  LEVELS  =’,I3, 
»  III) 

WRITE(6 ,203) 

203  FORMAT ( 1H  ,T5,’(  DATA  FORMAT  )»,//, T9, 'GAS#' ,T15 , 'WAVE  #',T24, 

»  'PRESSURE' ,T36, 'TEMP. ' , T46 , « PPM ’ , T58 , 'RANGE' ,T70, 'UGAS' ,T78, 

*  'TRAMSM. ' ,/) 
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DO  12  J=1 , JM 
WRITE (6 , 204 )  J 

204  FORMAT ( 1H0 , T5 , '  <  CUT’, 13,'  >’,/) 

T =0 . 

IT=0 

DO  13  K  =  1 , KM 

READ (5 , 103 )  KGAS , FREQ , RPRES , RTEMP , PPM , RANGE , RUGAS, TX 
103  FORMAT( 1 2 , F 1 0 . 3 »  E 1 1 . 4 , F9 . 3 , E 1 1 . 4 f  E 1 3 . 6 , E 1 1 .4.F7.4) 

RUGAS=RUGAS/CF 

WRITE (6 , 205 )  KGAS , FREQ , RPRES , RTEMP , PPM , RANGE , RUGAS , TX 

205  FORMAT (T 10, 12, FI 0.3, El  1 . 4 ,F9 . 3 , El  1 . 4 , El 3 . 6 , El  1 . 4 , F8 . 4) 

PRES,  TEMP  &  UGAS  ARE  CONVERTED  TO  THE  LOG  OF  THE  NORMALIZED 
VALUES.  IF  RPRES=Q  (INDICATES  NO  DATA), THEN  UGAS(M,J,K)  IS  SET 
AT  AN  IMPOSSIBLE  VALUE  ,  ALSO  RI(M,J,K)  IS  SET  TO  ZERO. 

IF ( R PRES . GT . 0 . )  GO  TO  14 
PRES(M,J,K)=0. 

TEMP (M, J , K) =0 . 

UGAS ( M , J , K) = 1 0 . 

RI (M, J , K) =0 . 

GO  TO. 15 

14  CONTINUE 

PRES (M,J,K)=AL0G10(R PRES/1  013. ) 

TEMP(M, J,K)=AL0G1 0(273. 15/RTEMP) 

UGAS(M,J,K)=ALOG 10 (RUGAS) 

15  CONTINUE 


T=T+TX 
IT=IT+1 
RI (M, J , K)= 1 .0 
13  CONTINUE 

TAU ( M, J) =T/ FLOAT ( IT ) 

12  CONTINUE 
GO  TO  1000 

END  OF  DATA  INPUT 

CONSTANTS  USED  IN  LATER  COMPUTATION  ARE  INITIALIZED 
FROM  2000  TO  3000. 

NCUT  =  MAXIMUM  //  OF  CUTS  USED  IN  COMPUTATION 
NDIM  =  DIMENSION  OF  THE  COEFFICIENT  MATRIX 

2000  CONTINUE 

IF(M.GT.O)  GO  TO  20 
WRITE (6 , 206 ) 

206  FORMAT(1H0,//,T1O, '$$$  NO  INPUT  DATA  $$$') 

STOP 

20  CONTINUE 
NC=M 

CSTD(1 )=0. 

NCUT  =NTC ( 1 ) 

NPTS=NTC ( 1 ) *NLV( 1 ) 

IF(NC.LE.I)  GO  TO  21 
DO  22  1=2, NC 
NCUT=MAX0(NCUT , NTC( I ) ) 
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N  PTS=NPTS+NTC (I)*NLV(I) 

22  CONTINUE 
21  CONTINUE 

FNCsFLOAT(NC) 
RIT=FLOAT(NPTS ) 

DO  23  J=1 ,NCUT 
TC=0 . 

DO  24  M= 1 , NC 
TC  =TC+TAU ( M , J ) 

24  CONTINUE 

TSTD( J )=TC/FNC 

23  CONTINUE 
NDIM=NC+1 +NCUT 


COMPUTATION  OF  THE  ABSORBER  PARAMETERS. 

THIS  LOOP  WILL  BE  REPEATED  MAXRPT  TIMES. 

FORMATION  OF  THE  NORMAL  EQUATION  AX  =  B  ,  A  IS  SYMMETRIC 

3000  CONTINUE 

DO  30.  1  =  1,19 
B ( I ) =0 . 

DO  31  J  =  1 ,19 
A ( I , J ) =0 . 

31  CONTINUE 
30  CONTINUE 

DO  1  M=1 , NC 
JM=NTC ( M) 

KM=NLV ( M) 

DO  2  J=1 , JM 
DO  3  K=1 ,KM 

IF(RI(M, J,K) .LT.0.5. )  GO  TO  3 

DO  4  IC=1 ,19 

V(IC)=0. 

4  CONTINUE 

V (NDIM-M ) = 1 . 

V(NDIM)  =  PRES (M , J , K) 

V (NDIM-1 )=TEMP(M,J,K) 

VV=-UGAS(M,J,K) 

V (NC UT+1 -J  )  =  1  . 

DO  5  11=1 , NDIM 
DO  6  IJ=1 , NDIM 

A(II,IJ)=V(II)*V(IJ)*RI(M,  J,K)  +  A ( II , I J ) 

6  CONTINUE 

B(II)=V(II)#VV*RI(M,J,K)  +  B ( II ) 

5  CONTINUE 
3  CONTINUE 
2  CONTINUE 
1  CONTINUE 

IF  J-TH  ROW  OF  "A"  IS  ZERO,  A(J,J)  IS  CHANGED  TO  -1 
WHICH  IS  DONE  IN  ORDER  TO  MAKE  "A"  NON-SINGULAR 
THIS  HAPPENS  WHEN  ALL  OF  THE  DATA  FOR  BAND  J-1  FAIL  TO  SATISFY 
THE  CRITERION  W  <  WCRIT.  THE  BAND  J-1  WILL  BE  IGNORED  IF  THIS 
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HAPPENS,  AND  THE  C-VALUE  FOR  BAND  J-1  WILL  BE  COMPUTED 
C  SEPARETELY. 

C 

ICONST ( 1 ) = 1 
IF ( NC . EQ . 1  )  GO  TO  40 
DO  41  M=2 , NC 
ICONST ( M) = 1 
I =NDIM-M 

IF( A( I , I ) . NE . 0. )  GO  TO  41 
A  ( I ,  I )  =  - 1 .0 
ICONST ( M) =0 

41  CONTINUE 
40  CONTINUE 

C 

C 

NC0L=0 

DO  42  J  =  1 , NDIM 
DO  43  1=1 , NDIM 
NC0L=NC0L+1 
X(NC0L)=A(I,J) 

43  CONTINUE 

42  CONTINUE 

PRINTING  OF  THE  HEADING  FOR  EACH  TRIAL  AND  THE  NORMAL  EQUATION 

IF ( LOOPCT . GT . 1 )  GO  TO  50 
WRITE (6 , 207  )  MAXR PT , LOOPCT 

207  FORMAT ( 1 H 1 , T20 ,' *#*  ABSORBER  PARAMETER  COMPUTATION  •«»',///, 

«  T 15 , 'NOTE :  THE  COMPUTATION  WILL  BE  REPEATED  MAXRPT  =',I2, 

*  '  TIMES. ',////, T10, 'TRIAL  # ' , I 1 , 5X , ' ( ALL  DATA  WERE  USED)') 

GO  TO  51 

50  CONTINUE 
WRITE (6 , 208  )  LOOPCT 

208  FORMATOH1  ,T10,  'TRIAL  1 1 , 5X ,’(  PARTIAL  DATA  WERE  USED  WITH’, 

«  '  CUT-OFF  CRITERION  :  W  <  2  )’) 

51  CONTINUE 

WRITE (6 , 209 )  NDIM, NDIM 

209  FORMAT(// , 1H0, ' <  NORMAL  EQUATION  :  AX  =  B  >',// ,T10 ,' .WHERE  THE', 

*  '  COEFFICIENT  MATRIX  A(»,I3,'  , ',13,'  )  AND  THE  CONSTANT  VECTOR1, 
«  '  B  ARE' ,//) 

IF(NDIM.LE. 17)  GO  TO  52 
WRITE (6 ,210)  NDIM 

210  FORMAT ( 1H  ,'***  WARNING  :  DIMENSION  OF  THE  MATRIX  IS  TOO  LARGE', 

«  '  (',13,'  )  TO  BE  PRINTED  IN  A  MATRIX  FORM  *»•',/) 

52  CONTINUE 
DO  53  1=1 , NDIM 

WRITE(6 ,211)  ( A( I , J ) , J  =  1 , NDIM)  ,B(I) 

211  FORMAT ( 1 H  , 1 8F7 - 3 ) 

53  CONTINUE 

****  MATRIX  INVERSION  SUBROUTINE  SIMQ  IN  SSP  IS  CALLED  ••*«••*§«•»•» 
CALL  SIMQ(X,B,NDIM,KS) 

PRINTING  OF  THE  SOLUTION  FOR  THE  NORMAL  EQUATION 
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IF(KS.EQ.I)  WRITE(6 ,212) 

212  FORMAT(1H0,T10, 'WARNING:  THE  COEFFICIENT  MATRIX  IS  SINGULAR.') 
AN=B(NDIM) 

AM=B(NDIM-1 ) 

IF(NC.LE.I)  GO  TO  54 
DO  55  M=2,NC 
CSTD (M)=B(NDIM-M) 

55  CONTINUE 
54  CONTINUE 

DO  56  J=1 ,NCUT 
PW( J)=-B(NCUT+1-J) 

56  CONTINUE 

WRITE (6, 21 3)  AN, AM, (CSTD(M) ,M=1 ,NC) 

213  FORMAT(// , 1H0, '  <  RESULTS  >',///, T7 N ', T1 7 M ', T27 ,' C 1 ', T37 ,' C2 • , 

*  T47, 'C3' ,T57, 'C4' ,T67, *C5  * ,T77, • C6  * , // , 2F 1 0 . 5 , 6F 10 . 3 ) 

WRITE (6,214)  (PW( I ) , 1  =  1 , NC  UT ) 

214  FORMAT (/, 1H0,T7, 'X*1 ’ , T 1 7 , ' X«2 • , T27 , ' X«3 ' , T37 , 'X«4 • , T47 , ' X«5 ’ , 

*  T57, ' X*6 ' ,T67, 'X*7 ' ,T77, '  X*8 ' ,T87, ' X *9  * ,  T97 , ' X*1 0 ' ,T107, 'X*1 1 ' , 

*  T1 17, 'X*12'  ,//, 12F10.3) 

NEL=NPTS-INT(RIT) 

WRITE(6 ,215)  NEL 

215  FORMAT (//, 1H0,T4, ' #  OF  ELIMINATED  POINTS  =',I5) 

CHECKING  OF  THE  CRITERION  (  W  <  WCRIT  )  AND  THE  COMPUTATION 
OF  C-VALUES  FOR  THE  IGNORED  BANDS. 

RI(M,J,K)  =  0  IF  W  IS  GREATER  THAN  OR  EQUAL  TO  WCRIT 
RI(M,J,K)  =  1  IF  W  IS  LESS  THAN  WCRIT 

RIT=0 . 

DO  60  M=1 ,NC 
JM=NTC ( M) 

KM=NLV(M) 

CAVG=0 . 0 
DO  61  J=1 , JM 
DO  62  K  =  1  ,  KM 

W=AN*PRES(M,J, K )+AM*TEMP ( M, J,K)+UGAS(M,J,K) 

IFCW.GE .WCRIT)  RI(M,J,K)=0. 

RIT=RIT+RI(M,J,K) 

CAVG=CAVG+ (PW( J )-W) 

62  CONTINUE 
61  CONTINUE 

IF( ICONST(M) . EQ. 1 )  GO  TO  60 
CSTD(M)=CAVG/FLOAT(JM«KM) 

WRITE(6 ,216)  M , M , M , CSTD( M) 

216  FORMAT ( /// , 1 H  ,T7,'»*  WARNING  »»',T25,'NO  DATA  FOR  BAND', 12, 

«  '  SATISFIES  THE  CRITERION  (  W  <  2  ).',//, T25 , 'THE  C '  ,11, 

*  ’  VALUE  IS  SEPARATELY  COMPUTED  BY  AVERAGING .',//, T30 ,' C ', I 1 , 

*  '  = ' ,F10. 3) 

50  CONTINUE 

COMPUTATIONS  OF  STANDARD  DEVIATIONS  IN  X 

NGDATA=0 
GTSD=0 . 

ICST=NC 
DO  70  M= 1 ,NC 
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JM=NTC ( M) 

KM=NLV(M) 

NDATA(M)=0 

TSD(M)=0. 

WRITE (6 ,201) 

WRITE (6, 202)  M,WN(M) ,NTC(M) ,NLV(M) 

WRITE(6 ,217)  AN , AM , M , CSTD(M) 

217  FORMAT(1H  ,T10,'N  = ’ ,F10.5,//,T10, 'M  = ' , F10 . 5, //,T10 , »C » , I 1 , 

*  '  = ' ,F 10 . 5 ) 

WRITE (6  218) 

218  FORMAT(//f 1H0.T7, 'RECOMPUTED  X-VALUES  AND  STANDARD  DEVIATIONS', 
»  '  IN  X-VALUES' ,/,1H0,T2, 'CUT' ,T11 , 'TAU' ,T20, 'X*' ,T30, 'XI ' ,T39, 

•  ' X2 ' , T48 , ' X3 ' ,T57, ' X4  * ,T66 , 'X5 » ,T75 , 'X6 ' ,T84 , *X7 ' , T93 , ’X8 ’ , 

*  T102, 'X9' , T 1 1 1 , 'X10' ,T121 , 'CUTWISE-SD' ,/) 

COMPUTATION  OF  THECUTWISE  STANDARD  DEVIATIONS  IN  X 

DO  71  J= 1 , JM 
DN=0. 

WW=0. 

DO  72  K  =  1  ,  KM 

P(K)=CSTD(M)+AN*PRES(M,J, K)+AM*TEMP( M ,J,K)+UGAS(M,J,K) 

WWW(J ,K)s(PW(J)-P(K) )*«2«RI(M, J,K) 

WW=WW+WWW( J , K ) 

DN=DN+RI (M , J , K) 

72  CONTINUE 
WW=SQRT ( WW/DN ) 

NDATA(M)=NDATA(M)+IFIX(DN) 

WRITE(6 ,219)  J,TAU(M,J),PW(J),(P(K),Ks1,KM) 

219  FORMAT ( 1H  ,15 ,?9. 3,F9. 4, IX, 10F9. *) 

WRITE(6 , 220 )  WW 

220  FORMAT (1H+,T121 , F 1 0 . 5 ) 

71  CONTINUE 

COMPUTATION  OF  THE  LEVELWISE  STANDARD  DEVIATIONS  IN  X 

DO  73  K  =  1 , KM 
WW=0. 

DN=0. 

DO  7M  J=1,JM 
WW=WW+WWW( J , K) 

DN=DN+RI (M, J , K) 

74  CONTINUE 

TSD(M) =TSD(M)+WW 
STANDV(K)=SQRT (WW/DN) 

73  CONTINUE 

WRITEC6, 221)  (STANDV(K ) ,K=1 , KM) 

221  FORMAT( 1H0,T4, 'LEVELWISE-SD  : • ,T26 , 10F9 . 5) 
NGDATA=NGDATA+NDATA(M) *ICONST(M) 
GTSD=GTSD+TSD(M)*FLOAT(ICONST(M) ) 

ICST=ICST-ICONST(M) 

TSD(M)=SQRT(TSD(M) /FLOAT (NDATA(M) ) ) 

WRITE (6 , 222 )  TSD(M) 

222  F0RMAT(//,1H0,T4, 'TOTAL  STANDARD  DEVIATION  FOR  THIS  BAND  :', 

•  F15.6) 

70  CONTINUE 
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PRINTOUT  OF  THE  SUMMARY. 

ALL  VITAL  INFORMATIONS  ARE  PRINTED  OUT  HERE. 


GTSD=SQRT (GTSD/ FLOAT (NGDATA ) ) 

WRITE (6 , 223 )  LOOPCT , AN , AM 

223  FORMAT ( 1H 1 , T15 , ’  ***  SUMMARY  OF  THE  ABSORBER  PARAMETER', 

*  •  COMPUTATION  FOR  TRIAL  #',12, '  •••',///, T20 , 

«  'PRESSURE  EXPONENT  N  = • ,F10 . 5 , // , T20 , 

*  'TEMPERATURE  EXPONENT  M  = • ,F 1 0 . 5 , /// ,T5 , ' CASE  #',3X, 

*  'WAVE  NUMBER'  ,5X,  'C -VALUE'  ,  5X ,  'TOTAL  if  OF  DATA',3X, 

*  ' CASEWISE  S.D.  IN  P') 

WRITE (6 , 224 )  (M,WN(M) ,CSTD(M) , NDATA( M) ,TSD(M) ,M=1 ,NC) 

224  FORMAT (1H0,T6,I3,6X,F9.2,5X,F8.3, 10X, 13, 12X.F12. 6) 

WRITE (6 , 225 )  NGDATA , NEL , GTSD 

225  FORMAT(//, 1H0,T15, 'GRAND  TOTAL  if  OF  DATA  = ’ , 15 , // , T15 , ’ #  OF', 

*  '  ELIMINATED  DATA  =', 15, // ,T15 , 'GLOBAL  STANDARD  DEVIATION  IN  P', 

*  '  =' ,F 12 . 6 , //) 

IF(ICST.LE.O)  GO  TO  75 
DO  76  M=1 , NC 

IF(ICONST(M) .EQ. 1 )  GO  TO  76 
WRITE (6 , 226  )  M 

226  FORMAT ( 1 H  ,T15,'NOTE:  THE  BAND', 13,’  IS  NOT  INCLUDED  IN  THE’, 

*  '  FINAL  STANDARD  DEVIATION') 

76  CONTINUE 

75  CONTINUE 

WRITE(6 , 227 )  LOOPCT, (TSTD(J) ,PW(J) ,J=1 ,NCUT) 

227  FORMAT (/// ,1H0,T15,'***  STANDARD  EMPIRICAL  TRANSMISSION’, 

*  '  FUNCTION  FOR  TRIAL  #', 12,’  *»*',///, T20 , 'TAU ', T35 ,' X« ',/ , 

*  (1H0,T17,F7.3,T30,F8.4)) 


IF ( RIT . GT . 0 . )  GO  TO  80 

IF  NO  INPUT  DATA  SATISFIES  THE  CRITERION,  THE  COMPUTATION  IS 
TERMINATED.  THE  MOST  RECENT  RESULTS  WILL  BE  USED  IN  THE  SEQUAL 

WRITE(6 , 228) 

228  FORMAT ( 1H1,//,T15, '$$$  NO  INPUT  DATA  SATISFIES  THE  CRITERION  OF' 

*  '  (  W  <  2  )  $$$' ,//,T15, '$$$  THE  COMPUTATION  FOR  THIS  STEP  IS' 

*  '  TERMINATED  $$$') 

GO  TO  4000 

80  CONTINUE 

L00PCT=L00PCT+1 

IF (LOOPCT . GT . MAXRPT )  GO  TO  4000 
GO  TO  3000 

4000  CONTINUE 

SUBROUTINE  COMPUTATIONS  FOLLOW 

IF ( INDXC 1 ) . LE.O)  GO  TO  90 
C 

CALL  NMBC 


noon  non 


90  CONTINUE 

IF(INDX(2 ) . LE. 0)  GO  TO  91 
CALL  INTPL1 

91  CONTINUE 

IF ( INDX( 3 ) . LE . 0 )  GO  TO  92 
CALL  INTPL2 

92  CONTINUE 
STOP 
END 

SUBROUTINE  NMBC 

COMPUTATION  OF  G'-VALUES  FOR  NON-MAJOR  BANDS 
DIMENSION  B( 15) ,CS(15) ,FS(15) 

COMMON  /PARM 1 /  TSTD ( 1 2 V ,PW ( 1 2 ) , WN (6 ) ,CSTD(6 ) , NCUT , NC , NAME (20 ) , 

»  AN,AM,CF,IC0NST(6) ,NEL 

WRITE(6 , 5 )  (NAME(I) ,1=1 ,20) 

5  FORMAT ( 1H 1 ,T15 ,20A4) 

WRITE(6 ,10) 

10  FORMAT(  1  HO ,  T15  ,  '  •••■  CALCULATION  OF  THE  SPECTRAL  PARAMETERS', 

»  'FOR  NON-MAJOR  BANDS  **»'///) 

DF=1 .E30 

11  CONTINUE 
NFREQ=0 

12  CONTINUE 
C=0. 

1=0 

15  CONTINUE 

READ  (5 , 20 )  KGAS ,  F'REQ ,  P ,  T ,  UGAS ,  TX 
20  FORMAT(I2,F10.3,E11 .4,F9.3,24X,E11 .4,F7.4) 

IF(KGAS.EQ.O)  GO  TO  25 
IF  (KGAS. LT. 0)  GO  TO  35 

THE  FOLLOWING  IF -STATEMENT  IS  INSERTED  TO  DETECT 
AND  TO  IGNORE  THE  INVALID  DATA  POINTS. 

IF (UGAS . GE . DF)  GO  TO  15 
1  =  1  +  1 
WX=FREQ 
UGAS=UGAS/CF 

C=C+(PW( I )-AN*AL0G10 (P/101 3 . )-AM*ALOG10 (273*  15/T)-AL0G10(UGAS) ) 
GO  TO  15 
25  G=C/FLOAT(I) 

NFREQ=NFREQ+1 
CS(NFREQ)=C 
FS(NFREQ)=WX 
DO  27  M= 1 , NC 

IF( ABS(WX-WN(M) ) .LE.O. 1 )  CS(NFREQ)=CSTD(M) 

27  CONTINUE 

IF(NFREQ. EQ. 10)  GO  TO  30 
GO  TO  12 
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30  CONTINUE 

WRITE (6, 31)  (FS(K) ,K=1 , NFREQ) 

31  F0RMAT( 1H0,2X, 'WAVE  NUMBER 2X , 10F 1 1 . 0) 

WRITE(6 , 32)  (CS (K ) , K=1 , NFREQ) 

32  FORMAT (1  HO ,5X, 'C  VALUES ', 2X , 10F 1 1 . 3//) 

GO  TO  11 

35  CONTINUE 

IF(NFREQ.EQ.O)  GO  TO  HO 
WRITE (6 ,31)  (FS(K) ,K=1 , NFREQ) 

WRITE(6 , 32)  (CS(K) ,K=1 , NFREQ) 

40  CONTINUE 
RETURN 
END 

SUBROUTINE  INTPL1 

COMPUTATION  OF  THE  STANDARD  PIECEWISE -ANALYTICAL  TRANSMISSION 
FUNCTION 

VERSION  1-1  «*  A3  CD  =  0  •* 

TAU  =  EXP(-10**(  A1 (I )+A2(I )*X  )) 

DIMENSION  SDCUT ( 1 5 ) ,  IC  UT  ( 1 5 ) , SDTC UT(15),ITCUT(15) 

COMMON  /PARM 1 /  TSTD( 1 2 ) , PW( 1 2 ) , WN (6 ) , CSTD(6 ) , NCUT , NC , NAME (20 ) , 

»  '*  AN ,  AM ,  CF ,  ICONST(6 )  ,  NEL 

COMMON  /PARM3/  A 1 ( 1 1 ) , A2 ( 1 1 ) , A3 ( 1 1 ) 

SSD=0. 

IT0TAL=0 

IM=NCUT-1 

JM=NCUT-2 


COMPUTATION  OF  THE  COEFFICIENTS  A 1 ( I ) ,  A2(I)  AND  A3 ( I ) 

CTX1=ALOG10(-ALOG(TSTD( 1  ))) 

DO  50  1=1 ,IM 
PDIF=PW(I)-PW(I+1 ) 

CTX2=AL0G10(-AL0G(TSTD(I+1 ))) 

A1(I)=(PW(I)*CTX2-PW(I+1 )»CTX1 )/PDIF 
A2(I )= (CTX1 -CTX2 ) /PDIF 
A  3 ( I )=0 . 

CTX1 =CTX2 

SDTCUT (I )  =  0 . 

50  ITCUT ( I  )  =  0 

THE  FIRST  AND  LAST  VALUES  OF  TSTD  AND  PW  ARE  CHANGED 
FOR  THE  TABLE  OUTPUT.  TRUE  VALUES  ARE  TEMPORARY  STORED 
IN  THE  RESERVE. 

TRES 1 =TSTD( 1 ) 

TRES2=TSTD(NCUT) 

PWRES 1 =PW( 1 ) 

PWRES2=PW(NCUT) 

TSTD(1 )=1 .0 
TSTD(NCUT)=0 . 0 
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PW(1 )=-1 .E70 
PW(NCUT )= 1 . E70 

PRINT  OUT  OF  THE  RESULTS 

WRITE (6 , 2 )  (NAME(I) ,1=1 ,20) 

2  FORMAT ( 1H 1 , T25 , 20 A 4 , // , T15 , '  PIECEWISE-ANALYTICAL  STANDARD*, 

*  'TRANSMISSION  FUNCTION ',//, T20 ,' TAU (X )  =', 

»  '  EXP(-10»*(  A 1  +  A2*X  )  )' ,///,T15, 'DATA: * ,T23, 'FROM  (TAU  ,', 

*  'X-VALUE)  TO  (  TAU  ,X-VALUE)  WITH  (  A1  ,  A2  ) ' ) 

WRITE (6, 3)  (TSTD(I) ,PW(I) ,TSTD(I+1 ) , PW( 1+1 ) , A1 ( I ) , A2(I ) , 1=1 , IM) 

3  FORMAT(1HO,T28,  '  ( '  ,F6. 3,  ’  ,',F7.3,’)  C,F6.3,’  , ' ,F7. 3, ' ) ' ,T74, 

WRITE  (6  \  4  )  (I,’WN(I),CSTD(I),I  =  1,NC) 

4  FORMAT(1HO,//,T15, 'ABSORPTION  BANDS: ',T40, 

*  •#  WAVENUMBER  C-VALUE ' /( 1H0, T39 , 12 , 5X ,F7 . 1 ,F1 1 . 5) ) 

CHECK  IF  ANY  DATA  IS  AVAILABLE  FOR  S.D.  COMPUTATION 
DATA  FORMAT  IS  THE  SAME  AS  THAT  FOR  MAIN  PROGRAMME 
ONE  CONTROL  CARD  IS  READ-IN  FIRST  FOR  BRANCHING 
IFQ  >  0  DATA  SET  FOLLOWS,  READ-IN  DATA 

IFQ  =  0  END  OF  DATA,  GO  TO  THE  FINAL  PRINTING 

READ (5 ,11,END=42)  FQ, IFQ 

11  FORMAT (5X,F10.3,T41,I4) 

IF(IFQ.GT.O)  GO  TO  18 

42  WRITE(6 , 41 ) 

41  FORMAT (///,1HO,T5,'$$$  NO  DATA  FOR  STANDARD  DEVIATION  COMPUTATION 

*  $$$') 

GO  TO  40 

8  READ(5, 11 ,END=30)  FQ, IFQ 
IF(IFQ.LE.O)  GO  TO  30 
18  CONTINUE 
ST=0 . 

DO  51  1  =  1  ,IM 
SDCUT (I )=0. 

51  ICUT ( I ) =0 
CLOG= 1 00 . 

DO  52  1=1 ,NC 

52  IF ( ABS ( FQ-WN ( I ) ) . LT . 1 . )  CLOG=CSTD(I) 

IF ( CLOG . LT. 99 . )  GO  TO  13 

THE  READ-IN  WAVENUMBER  DOES  NOT  MATCH  THE  MAJOR  BAND 
WAVENUMBER  (WN(I)).THE  DATA  IN  THIS  BAND  ARE  IGNORED. 

WRITE (6 , 12 )  FQ 

12  FORMAT ( 1  HO , T10 , '**  ERROR  IN  WAVENUMBER  »»', 

•  '  (READ-IN  WAVENUMBER  =',F10.5,’  )') 

DO  61  IDUM= 1 , IFQ 
READ(5 , 60 )  DUMMY 

60  FORMAT(FI.O) 

61  CONTINUE 
GO  TO  8 

VALID  DATA  INPUTS , READ-IN  OF  THE  DATA  AND  STANDARD 
DEVIATION  COMPUTATION  ARE  PERFORMED  SIMULTANEOUSLY. 
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13  CONTINUE 
WRITE(6 ,17)  FQ 

17  F0RMAK1H1  ,T15,  ’  (  WAVE  NUMBER  =',F8.1,'  )',///,  6X, ' WAVEN .', 3X , 

*  'PRESS. ' ,4X, 'TEMP. ' ,7X, 'U • ,8X, 'TRANSM.  -  T(COMP)  =  DIFF',6X, 

*  ' DIFF**2  ' ,  4X , ' X -VALUE ' , / ) 

DO  9  M= 1 , NDATA 

READ(5, 10)  KGAS.FQ, PRES, TEMP, UG,TX 
10  FORMAT d2, F 10 . 3 ,  E 1 1  •  4 ,F9 . 3 , 24X,  El 1.4.F7.4) 

UG=UG/CF 

X=CLOG+AN*ALOG 10 (PRES/1013. )+AM*ALOG 10 (273 • 1 5/TEMP )+ALOG 1 0 (UG ) 

DO  14  J  =  1  ,  JM 

IF ( X . LE . PW( J+1 ))  GO  TO  15 

14  CONTINUE 
J  =  IM 

15  TC=EXP(-10**(A1 (J)+A2(J)*X)  ) 

D=TX-TC 
SD=D*D 
ST =ST+SD 

SDCUT( J )=SDCUT ( J )+SD 
ICUT(J)=ICUT(J)+1 

WRITE (6, 16)  FQ, PRES, TEMP, UG , TX , TC , D, SD ,  X 

16  FORMAT ( 1 H  ,3X,F8.1,F10.2,F9.2,E13.4,F9.4,F12.4,F13.6,E12.3,F9.3> 

9  CONTINUE 

END  OF  DATA  READ-IN  FOR  THIS  BAND. 

TOTAL  STANDARD  DEVIATIONS  ARE  COMPUTED  AND  PRINTED. 

20  SSD=SSD+ST 
ITOTAL=ITOTAL+IFQ 
ST=SQRT(ST/FLOAT(IFQ) ) 

DO  21  1=1 ,IM 

SDTCUT (I )=SDTCUT(I )+SDCUT(I ) 

ITCUT ( I ) =ITCUT(I )+ICUT (I ) 

21  SDCUT ( I ) =SQRT (SDCUT ( I ) /FLOAT (ICUT(I))) 

WRITE (6, 22)  (I , TSTD(I ) ,TSTD(I+1 ) , ICUT (I ) , SDCUT (I ) , 1=1 ,IM) 

22  FORMAT  ( 1 H0-, ///,  T 1 0  ,  'CUTWISE  STANDARD  DEVIATION  ',//,  T15  ,  *  T20  , 

*  '(  FROM,  TO  )  ' , T40 , ' #  OF  DATA ', T53 CUTWISE  SD • , // , (T 1 4 , 12 , 

*  T20, ’ ( ’ ,F5.2, ’ , ’ ,F5.2, ’ ) » , T43 , 14, T52 , FI 0 . 6 , /) ) 

WRITE (6 ,23 )  IFQ, ST 

23  FORMAT(1H0,T10, 'TOTAL  //  OF  DATA  FOR  THIS  BAND  =',I5,9X, 

«  '  STANDARD  DEVIATION  =',F12.6,//) 

GO  TO  8 

END  OF  THE  STANDARD  DEVIATION  COMPUTATION  FOR  ALL  DATA. 

GRAND  TOTAL  STANDARD  DEVIATION  IS  COMPUTED  AND  PRINTED  OUT 
TOGETHER  WITH  VITAL  INFORMATIONS. 

30  SSD=SQRT (SS D/FLO AT ( ITOTAL) ) 

DO  31  1=1 ,IM 

31  SDTCUT  (I  )  =  SQRT  (SDTCUT (I )  /  FLOAT dTC UT d ) ) ) 

WRITE (6, 32)  (I , TSTD(I ) , PW( I ) ,TSTD(I+1 ),PW(I+1),A1(I),A2(I), 

*  ITCUT(I) ,SDTCUT(I) ,1=1 ,IM) 

32  FORMAT (1H1 ,T20, ’ •**  PIECEWISE -ANALYTICAL  STANDARD  TRANSMISSION', 
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•  ’  FUNCTION  ,///,T10, 'TOTAL  CUTWISE  STANDARD  DEVIATION’,//, 

•  T15, 'CURVE  # ' , 3X  , 'FROM  (  TAU  .X-VALUE)  TO  (  TAU  .X-VALUE)', 

•  ’  WITH  (  A 1  ,  A2  )  '  ,3X,  '#  OF  DATA’ , ^X , 'CUTWISE  SD',//, 

•  (T18,I2,T30, •( ’ ,F6.3, ’  ,',F7.3,')  (',F6.3,'  , ' ,F7 . 3 , ’ ) ' , T76 , 

«  '(’ ,F7.4,’ ,F7.4, ') ' ,7X,I3,5X,F10.6,/)) 

WRITE (6 , 33 )  ITOTAL , SSD 

33  FORMATOHO.TIO, 'GLOBAL  RESULTS ',//, T15 , 'TOTAL  NUMBER  OF  DATA’, 

•  '  USED' ,15,//, T15, 'GLOBAL  STANDARD  DEVIATION' ,F12. 6) 

40  CONTINUE 

END  OF  ALL  COMPUTATION. 

RESERVED  TRUE  VALUES  OF  THE  FIRST  AND  LAST  TSTD 
AND  PW  ARE  RETURNED. 

TSTD( 1 )=TRES1 
TSTD(NCUT)=TRES2 
PW( 1 ) =PWR£S 1 
PW(NCUT)=PWRES2 

CALL  SDTAU 

RETURN 

END 

SUBROUTINE  INTPL2 

COMPUTATION  OF  THE  PIECEWISE-ANALYTICAL  STANDARD  TRANSMISSION 
FUNCTION 

VERSION  2  -  1 

TAU=EXP ( - 1 0** ( A 1 +A2*X+A3*X**2 ) ) 

COMMON  /PARM 1 /  TSTD( 12),PW(12),WN(6), CSTD(6 ) , NCUT, NC , NAME (20 ) . 

»  AN , AM , CF , ICONST(6 ) ,NEL 

COMMON  /PARM3/  A 1 ( 1 1 ) , A2 ( 11 ) , A3 ( 1 1  ) 

DIMENSION  TOO)  ,TD(6,74)  ,  PD  (6, 74)  ,  JI(6, 10)  , 

•  SDK (7 ) , SDE (7 , 9 ) , SUME(2 , 9 ) , DE (6 , 9 ) 

NWC=0 

K=0 

READ(5 , 2 , END=80 )  FREQ, MAXDAT 

ifcmaxdat;gt.o)  go  to  3 

80  CONTINUE 
WRITE(6 , 99 ) 

99  FORMAT (///, 1H0.T5, '$$$  NO  DATA  FOR  STANDARD  DEVIATION  COMPUTATION 

•  $$$') 

GO  TO  77 

1  CONTINUE 

READ(5,2,END=21 )  FREQ, MAXDAT 

2  F0RMAT(5X ,F 10. 5 , T41 ,14) 

IF(MAXDAT.EQ.O)  GO  TO  21 

3  CLOG= 1 . E  10 
DO  5  Ls1,NC 

IF( ABSCFREQ-WN(L) ) .LE.  .01)  CLOG=CSTD(L) 

5  CONTINUE 

IF ( ABS ( CLOG- 1 . E 1 0 ) . GE . . 0 1 )  GO  TO  9 
WRITE(6 ,100)  FREQ 
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100  FORMAT ( '1 ' ,///////, '  ERROR  IN  INPUT  DATA  ;  WAVE  NUMBER  ’,F10.3, 

*  ’  NOT  USED  IN  COMPUTATION  OF  CONSTANTS.') 

DO  6  J  =  1 , MAXDAT 

READ(5 , 101 )  KGAS,  FREQ, PRES, TEMP, PPM, RANGE, UGAS.TX 
6  CONTINUE 
GO  TO  1 
9  CONTINUE 
K=K+1 
NWC=NWC+1 
JI(K, 1 ) =0 
JI  (K  ,  NCUT )=MAXDAT 
J =2 

DO  20  1=1, MAXDAT 

READ(5, 101 )  KGAS, FREQ, PRES, TEMP, PPM, RANGE, UGAS,TX 

101  F0RMAT(T2,F10.3,E11 . 4 , F9 . 3 , El  1 .4, El  3. 6, Ell .4.F7.4) 

TD(K , I )=TX 

PD(K,I)=  AN  *ALOG10(PRES/1013.)+  AM  «ALOG 1 0 (273 . 1 5/TEMP )+ALOG 1 0 

*  ( UGAS/CF ) +CLOG 

IF(TD(K,I) .GE.TSTD(J) )  GO  TO  20 
IF(J.EQ.NCUT)  GO  TO  20 
JI  (K,  J)  =  I-1 
J  =  J  +  1 

20  CONTINUE 
GO  TO  1 

21  CONTINUE 
IF(NWC.LE.O)  RETURN 
DO  30  J= 1 , NCUT 

T ( J  )  =  ALOG 1 0 (-ALOG (TSTD ( J ) ) ) 

30  CONTINUE 
SUMT=0 . 

DT=0 . 0 
NCC=NCUT-1 
DO  45  1=1 . NCC 
SA=0  . 

TA=0  . 

UA=0 . 

DO  41  K= 1 , NWC 
SUME (K , I ) =0 . 0 
M  =  JI (K  ,  I )  +  1 
N  = JI (K , 1  +  1 ) 

DO  40  J  =M , N 

TC=ALOG 10 (-ALOG (TD (K , J ) ) ) 

S A  =SA+  TC  «PW(I)*PW(I  +  1)«(PD(K,J)-PW(I))*(PD(K,J)-PW(I  +  1)) 
TA=TA+((PD(K, J)-PW(I))*(PD(K, J)-PW(I+1 ))«( (PD(K, J>**2 >• (PWCI+1 ) 

*  *T(I)-PW(I)*T(I+1))+PD(K,J)*((PW(I)**2)*T(I+1 )- (PW ( 1+1 ) *»2 ) «T (I ) 

*  )  )  )/(PW(I)-PW(I  +  1 )) 

U A  =  U A  +  ( ( PD (K , J )-PW( I ) ) * (PD (K , J )-PW( 1  +  1 ) ) ) **2 
4  CONTINUE 
:  5NTINUE 

A ’  I )= (SA-TA)/UA 

A ‘ f ; ' s ( T ( I ) -T ( I ♦ 1 ) ) /( PW ( I + 1 )*(PW( I )-PW( 1+1 )))-(T(I)-A1(I))/(PW(I)* 

*  ?  d ! :  *  i ) ) 

»  •  -S (T(I )-T(I+1 ))/(PW(I)-PW(I+1 ))-A3(I)*(PW(I)+PW(I+1)) 
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M= JI (K , I )+ 1 
N  =  JI (K , 1+1 ) 

DE (K , I )=FLOAT ( 1 +N-M) 

DO  43  J=M,N 

SUME(K,I)  =  SUME(K,I)+(TD(K,J)-EXP(-10.**(A3'(I)*PD(K,J)*PD(K,J)  + 

«  A2(I)*PD(K,J)+A1(I))))**2 

43  CONTINUE 

SDE(K,I)=SQRT(SUME(K,I)/DE(K,I)) 

SUMI=SUMI+SUME(K , I )* FLOAT (ICONST(K)) 

DI =DI+DE (K ,  I ) * FLOAT (ICONST(K)) 

44  CONTINUE 
SUMT=SUMT+SUMI 
DT=DT+DI 

SDE (NWC+1 ,I)=SQRT(SUMI/DI) 

45  CONTINUE 

DO  51  K= 1 , NWC 
SUMK=0.0 
DK=0 . 0 

DO  50  1=1 ,NCC 
SUMK=SUMK+SUME(K,I) 

DK=DK+DE(K, I) 

50  CONTINUE 

SDK(K ) =SQRT (SUMK/DK ) 

51  CONTINUE 

SDK(NWC+1 )=SQRT(SUMT/DT) 

DUM 1 =PW ( 1 ) 

DUM2=PW ( NCUT ) 

DUM3=TSTD(1 ) 

DUM4=TSTD(NCUT) 

PW( 1 >=-1000000. 

PW (NCUT >=100 00 00. 

TSTDC1 )=1 .0 
TSTD(NC  UT )  =  0 . 0 
DO  60  K= 1 , NWC 

WRITE (6 , 102) (NAME(J) , J=1 ,20) ,WN(K) ,NCC, (I,TSTD(I) ,TSTD(I+1 ) ,PW(I) , 

*  PW(I+1),A1(I),A2(I),A3(I), SDE (K,I) ,1=1 ,NCC) 

102  FORMAT  ('  1  35X ,  ' RENDITION  OF  EMPIRICAL  TRANSMITTANCE  FUNCTION 

•FOR  :  ' // , 20A4 , // , 40X , 'WAVE  NUMBER  :  '  ,F 15 . 4 , //// , 20X , 'THE  TRANSMI 
•SSION  CURVE  IS  DIVIDED  INTO', 13,'  SEPARATE  CURVES .', /20X ,’ EACH  CU 
«RVE  IS  EXPRESSED  BY  A  FUNCTION  OF  THE  FORM  "  TAU  =  EXP(-1 0««(A3*P 
#*P+A2*P+A 1 ) )  " . ' / , 20X ,  'THE  FUNCTION  COEFFI 

•CIENTS  AND  RESULTING  STANDARD  DEVIATION  FOR  EACH  CURVE  ARE  AS  FOL 
•LOWS: ' ,//// 22X, 'TAU' , 20X , ' P ' , 24X , ' A 1 ' , 1 3X , ' A2 * , 1 3X , ' A3 ' , 1 1 X , 'STAND 
•ARD  DEVIATION',///,  (5X, 'CURVE  # ' , 13, 3X, • ( ' ,F4. 2, ' - • ,F4. 2, ' ) • ,5X, 
*' ( ' ,F9.5, '-' ,F9.5, ' ) ' ,5X,3F15.6,  6X,F15.6  /) ) 

WRITE (6 ,401)  SDK (K) 

401  FORMAT(1X,//,87X, 'TOTAL  STANDARD  DEVIATION ', F 1 5 . 6 ) 

60  CONTINUE 
K=K  + 1 

WRITE ( 6 , 1 04 )  (NAME(J) , J  =  1 ,20)  ,NCC,  ( I , TSTD( I ) , TSTD( 1  +  1 ) , PW( I ) , 

•  PW ( 1  +  1 ) , A 1 (I ) ,A2(I ) , A  3 ( I ) i SDE (K,I),I  =  1,NCC) 

104  FORMAT  (» 1 »,/, 35X ,» RENDITION  OF  EMPIRICAL  TRANSMITTANCE  FUNCTION 
•FOR  :  ’//,20A4,//,40X, 'TOTAL  PROFILE  AVERAGED  OVER  ALL  WAVE  NUMBER 
•S'  ,////, 20X, 'THE  TRANSMI 

•SSION  CURVE  IS  DIVIDED  INTO', 13,'  SEPARATE  CURVES. ' ,/20X , 'EACH  CU 
•RVE  IS  EXPRESSED  BY  A  FUNCTION  OF  THE  FORM  "  TAU  =  EXP(-10«« (A3*P 


■ 
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#*P+A2*P+A 1 ) )  " . ' / , 20X ,  'THE  FUNCTION  COEFFI 

•CIENTS  AND  RESULTING  STANDARD  DEVIATION  FOR  EACH  CURVE  ARE  AS  FOL 
•LOWS: ' ,////22X, ’TAU ’ , 20X , ’ P ' , 24X, • A 1 ' , 1 3X, ' A2 ’ , 1 3X, • A3 ' , 1 1X, 'STAND 
•ARD  DEVIATION’,///,  (5X, 'CURVE  # ’ , 13, 3X, ’ ( ' ,F4 . 2, ' - ' ,F4. 2, ' ) ' ,5X , 
•' ( ' ,F9.5, '-' ,F9.5, ’ ) ' ,5X,3F15.6,  6X.F15.6  /)) 

WRITE(6 , 402  )  SDK (K ) 

402  F0RMAT(1X,//,81X, 'GRAND  TOTAL  STANDARD  DEVIATION ' ,F15 . 6 ) 

WRITE  (7,201 )  (AHI)  ,A2(I)  ,A3(I),I  =  1  ,NCC) 

201  FORMAT (3F 1 0 . 6 ) 

IDT=IFIX (DT ) 

WRITE (6, 225) IDT, NEL, SDK (K) 

225  F0RMAT(//,  1H0,T15, 'GRAND  TOTAL  OF  DATA  a  •  ,  15,//,  T15  , '#  OF', 

*  ’  ELIMINATED  DATA  =', 15, //, T15 , 'GLOBAL  STANDARD  DEVIATION  IN', 

*  ’  TAU  =' ,F 12 . 6 , // ) 

DO  76  M= 1 , NWC 

IF ( ICONST (M) . EQ . 1 )  GO  TO  76 
WRITE (6 , 226  )  M 

226  FORMAT ( 1 H  ,T15,'NOTE:  THE  BAND', 13,’  IS  NOT  INCLUDED  IN  THE’, 

*  '  FINAL  STANDARD  DEVIATION') 

76  CONTINUE 
PW( 1 ) =DUM 1 
PW ( NC  UT ) =DUM2 
TSTD(1 )=DUM3 
TSTD(NCUT)=DUM4 

CALL  SDTAU 

77  CONTINUE 
RETURN 
END 

SUBROUTINE  SDTAU 

COMPUTATIONS  OF  STANDARD  DEVIATIONS  IN  TAU  USING  THE  ORIGINAL 
DATA  USED  IN  MAIN 

DIMENSION  NDATA( 6 ) ,TSD(6 ) ,WWW( 12,10), ST AND V (12),P(10),T(10) 

COMMON  /PARM 1 /  TSTD( 12) , PW( 1 2 ) ,WN(6) ,CSTD(6) ,NCUT,NC, NAME (20) , 

*  AN, AM,CF, ICONST (6) ,NEL 
COMMON  /PARM2/  PRES ( 6 , 1 2 , 1 0 ) , TEMP ( 6 , 1 2 , 1 0 ) , UGAS (6 , 1 2 , 1 0) , 

*  TAU(6 , 12 ) , NTC (6 ) , NLV(6 ) 

COMMON  /PARM3/  A1 ( 1 1 ) , A2( 1 1 ) , A3( 1 1 > 

C 

NGDATA=0 
GTSD=0 . 

ICST=NC 
DO  70  M= 1 , NC 
JM=NTC ( M) 

KM=NLV(M) 

NDATA ( M) = JM*KM 

NGDATA =NGDATA+N DATA ( M) * ICONST (M) 

TSD(M) =0 . 

WRITE(6 ,214) 

214  FORMAT( ’ 1 ’,/////, 45X, 'RECOMPUTATION  OF  TAU',/////) 

IF(M.GT.I)  GO  TO  77 
WRITE(6 ,215) 


o  o  o  o  o  o 
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215  FORMAT(20X, 'A  TAU  VALUE  ,  T  ,  IS  RECOMPUTED  FOR  THE  ORIGIONAL  DATA 
»  USING  THE  PIECEWISE-ANALITICAL  TRANSMISSION  FUNCTION //20X, 

*  'STANDARD  DEVIATIONS  BETWEEN  THE  ACTUAL  TAU  AND  THE  RECOMPUTED', 

*  '  TAU  VALUES  ARE  COMPUTED.*////) 

77  CONTINUE 

WRITE(6 ,202)  M,WN(M) ,NTC(M) ,NLV(M) 

202  F0RMAT(1H0,T15, '•••  CASE', 13,'  (WAVE  NUMBER  =' ,F10. 3 ,' ) 

*  ///,T20, 'TOTAL  #  OF  CUTS  s ', 13 ,// ,T20 , 'TOTAL  #  OF  LEVELS  =',I3, 

*  ///) 

WRITE(6 ,216)  AN , AM ,M ,CSTD(M) 

216  F0RMAT(10X, 'N  =  ' , F 1 0 . 5 , // 1 0X ,  ' M  = ' ,F 10 . 5 , // ,  10X , ' C ' , 1 1 , '  =', 

*  F10.5,////) 

WRITE(6 ,217) 

217  FORMAT(//,1HO,T7, 'RECOMPUTED  TAU  AND  STANDARD  DEVIATIONS', 

*  '  IN  TAU  ',/,1H0,T2,'CUT',T11,'TAU',T20,'X«',T30,'X1’,T39, 

*  'X2' ,T48, *X3' ,T57, 'X4* ,T66, 'X5' ,T75, 'X6',T84, 'X7' ,T93, ’X8’ , 

*  T102, 'X9' , T 1 1 1 , ' X 1 0 ’ , T 1 2 1 , 'CUTWISE-SD' ,/) 

COMPUTATION  OF  THE  CUTWISE  STANDARD  DEVIATIONS  IN  X 

DO  71  J  =  1 , JM 
WW=0. 

DO  72  K= 1 , KM 

P(K)=CSTD(M)+AN*PRES(M,J, K )+AM*TEMP( M , J , K)+UGAS( M , J , K ) 

IM=JM-1 
DO  75  1=1, IM 

IF( PW( 1+1 ) .GT.P(K) )  GO  TO  76 

75  CONTINUE 
I  =  IM 

76  CONTINUE 

T(K)=EXP(-10*«(A3(I)*P(K)«P(K)+A2(I)»P(K)+A1 (I))) 

WWW( J ,K)=(TAU(M,J)-T(K))**2 
WW=WW+WWW( J , K) 

72  CONTINUE 

WW=SQRT (WW/ FLOAT (KM ) ) 

WRITE (6, 218)  J,PW(J) ,TAU(M, J) ,(T(K) ,K=1 ,KM) 

218  FORMAT ( 1 H  , 15 , F9 . 4 , F9 . 4 , 1 X , 10F9 . 4 ) 

WRITE (6 ,219)  WW 

219  FORMAT (1H+,T121 , F 1 0 . 5 ) 

71  CONTINUE  . 

COMPUTATION  OF  THE  LEVELWISE  STANDARD  DEVIATIONS  IN  X 

DO  73  K  =  1 , KM 
WW=0. 

DO  74  J=1 , JM 
WW=WW+WWW(J ,K) 

74  CONTINUE 

TSD(M) =TSD(M)+WW 

STANDV(K )=SQRT (WW/ FLOAT ( JM ) ) 

73  CONTINUE 

WRITE(6 ,220)  (STANDV(K) ,K=1 ,KM) 

220  FORMAT( 1H0, T4, 'LEVELWISE-SD  : ' ,T26, 10F9. 5) 

GTSD=GTSD+TSD(M) *FLOAT (ICONST(M) ) 

IC ST=ICST-I CONST (M) 

TSD(M)=SQRT(TSD(M) /FLOAT (NDATA(M) ) ) 
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WRITE(6 ,221 )  TSD(M) 

221  F0RMAT(1H0,//,T15, 'TOTAL  STANDARD  DEVIATION  FOR  THIS  CASE  :', 

•  F15.6) 

70  CONTINUE 

GTSD=SQRT(GTSD/FLOAT(NGDATA  ) ) 

WRITEC6  223)  AN  AM 

22 3  FORMAT ( '1 ' ,T15, '***  SUMMARY  OF  THE  TRANSMITTANCE  RECOMPUTATION  * 
«  **’ ,///,T20, 'PRESSURE  EXPONENT  N  = * ,F 10 . 5 , // , T20 , 

*  'TEMPERATURE  EXPONENT  M  = ' ,F 1 0 . 5 , /// , T5 , ' CASE  #',3X, 

«  'WAVE  NUMBER' ,5X, 'C-VALUE' ,5X, 'TOTAL  #  OF  DATA',3X, 

#  'CASEWISE  S.D.  IN  TAU') 

WRITE (6, 224 ) ' (M,WN(M) , CSTD( M) , NDATA( M) ,TSD(M) ,M=1 , NC) 

224  FORMAT(1HO,T6,I3,6X,F9.2,5X,F8.3, 10X,I3, 12X.F12.6) 

WRITE (6 , 225  )  NGDATA , NEL , GTSD 

225  FORMAT(//,1HO,T15, 'GRAND  TOTAL  #  OF  DATA  = ' , 15 , // , T1 5 , ' #  OF', 

«  '  ELIMINATED  DATA  =', 15 ,//, T15 , 'GLOBAL  STANDARD  DEVIATION  IN’, 

*  '  TAU  =’ ,F 12 . 6 , //) 

IF(ICST.LE.O)  RETURN 
DO  78  M= 1 , NC 

IF(ICONST(M) .EQ. 1 )  GO  TO  78 
WRITE (6 , 226  )  M 

226  FORMAT ( 1H  ,T15,'N0TE:  THE  BAND’, 13,'  IS  NOT  INCLUDED  IN  THE', 

»  '  FINAL  STANDARD  DEVIATION') 

78  CONTINUE 
RETURN 
END 


ooooooooooaoooooooooonoooooooaoanoooooooo 
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COMPUTER  CODE  SIMMIN  ••»••• 

VERSION  (  6  -  3  )  TRACE  GASSES 

COMPUTATION  OF  ABSORBER  PARAMETERS  AND  ANALYTICAL  STANDARD 
TRANSMISSION  FUNCTION 


THIS  CODE  USES  THE  SUBROUTINE  FMCG  IN  SSP  LIBRARY 
THIS  CODE  CONSISTS  OF 

MAIN:  DATA  READ-IN  AND  CONTROL  OF  COMPUTATION 

FUNCT:  COMPUTATION  OF  THE  COST  AND  ITS  DERIVATIVES 

TITLE:  PRINTING  OF  HEADINGS  AND  INITIAL  CONDITIONS 

PRINT1 :  PRINTOUT  OF  RESULTS  AND  COMPUTATION/PRINTING  OF  S.D.S 

NMBC :  COMPUTATION  OF  NON-MAJOR  BANDS'  C-VALUES 


DATA  SET-UP 

1.  INITIAL  GUESSES  X(I)  (9  CARDS  WITH  T12.F10.7) 

X(  1 )  =A  1 ,  X  ( 2  )  =A2 ,  X(3)=A'3,  X(4  )sN  ,  X(5)=M,  X(5+I  )=LOG(C(  I ) ) 
(NEED  DUMMY  INPUTS  FOR  PROBLEMS  WITH  DIMENSION  <  9) 

2.  SIGNAL  VARIABLES  S(I)  (9  CARDS  WITH  T12.F10.7) 

S(I)  =0  ...  X(I)  IS  KEPT  CONSTANT 
S ( I )  =  1  ...  X(I)  IS  VARIED 

(NEED  DUMMY  INPUTS  FOR  PROBLEMS  WITH  DIMENSION  <  9) 

3.  COMMENT  CARD  (20A4)  FOR  TITLE  AND  ABSORBER  TYPE  ETC. 

4.  DATA  SETS  (MAX.  4  SETS)  -  ONE  FOR  EACH  ABSORPTION  BAND 

EACH  SET  CONSISTS  OF 

IST(CONTROL)  CARD:  WAVENUMBER,  #  OF  DATA. AND  COMMENTS 

(SEE  FORMAT  101 ) 

DATA  CARDS:  P,  T,  U,  TAU  ETC. 

(SEE  FORMAT  102) 

(TOTAL  #  OF  DATA  SHOULD  NOT  EXCEED  900) 

5.  BLANK  CARD  -•  FOR  THE  TERMINATION  OF  DATA  INPUT  FOR  MAIN 

6.  DATA  SETS  FOR  NMBC  -  ONE  FOR  EACH  ABSORPTION  BAND 

EACH  SET  CONSISTS  OF 

DATA  CARDS:  SAME  AS  MAIN 
FINAL  CARD:  BLANK 

TERMINATION  A  CONTROL  CARD  WITH  -1  IN  FIRST  TWO  COLUMNS 
THIS  COMES  AFTER  THE  FINAL  BLANK  CARD 
(IF  NO  DATA  BUT  A  BLANK  CARD  IS  SUPPLIED,  NMBC  IS  SKIPPED) 


DIMENSION  X(9) ,G(9) ,Y(9) ,H(72) ,WN(4) 

COMMON  /PARM1 /  NC , ND ( 5 )  ,  RW ( 4  ) 

COMMON  /PARM2/  IC , PLOG ( 900 ) , TLOG ( 900 ) , ULOG ( 900 ) , TAU (900 ) , S(9 ) 
COMMON  /PARM3/  P ( 900 ) , T ( 900 ) , U (900 ) , L(20 ) 

COMMON  /PARM4/  PO , TO , NDIM , ID (5 , 9 ) 

EXTERNAL  FUNCT 
C 

C  CONSTANTS 

P0=1 .013E+03 
T0=273- 15 
CF=2 . 69E+1 9 
N  =  9 
V  =  0. 

IC=0 

MAXNC=4 


ooo  o  non  oooooo  oo  oo 


DATA  INPUT 

READ(5 ,100)  (X(I),I=1,N) 

READ(5 ,100)  (S (I ) , 1= 1 , N) 

100  FORMAT  (T12.F10.7) 

COMMENT  CARD  (THIS  INCLUDES  THE  ABSORBER  TYPE) 

READC5 , 500)  (L ( I ) , 1= 1 , 20 ) 

500  FORMAT(20A4 ) 

NC  =  #  OF  MAJOR  ABSORPTION  BANDS 

ND(1)=0,  ND(2)=N1,  ND(3 )=N 1+N2,  ND(4 )=N 1+N2+N3 ,  ... 

WHERE  N 1 ,  N2 ,  N3,  ...  ARE  #S  OF  DATA  IN  BANDS  1,  2,  3,  ... 
NDCNC  +  1  )  =  TOTAL  if  OF  DATA 

NC=0 
ND ( 1  )  =  0 

DO  10  M  =  1  , MAXNC 

READC5, 101)  WN(M) , IX , ( ID (NC+1 , I ) ,1=1 ,9) 

101  FORMAT (5X,F10.3,T41,I4,9A4) 

IF(IX.LE.O)  GO  TO  11 
NC=NC+1 

IM=ND (NC ) + 1 
IN=ND(NC)+IX 
ND(NC+1 )=IN 
DO  12  IaIM.IN 

READ(5 ,102)  P(I),T(I),U(I) ,TAU(I ) 

102  FORMAT  ( 1 2X, El  1 . 4 ,F9. 3 ,24X , El  1 . 4,F7 . 4) 

U(I)=U(I)/CF 

DATA  ARE  CONVERTED  TO  THE  LOG  OF  THE  NORMALIZED  VALUES 

PLOG ( I ) = ALOG 10(P(I)/P0) 

TLOG ( I ) = ALOG 10(T0/T(I)) 

ULOG (I )  =  ALOG 10 (U ( I ) ) 

12  CONTINUE 

10  CONTINUE 

11  CONTINUE  ' 

END  OF  DATA  INPUT 

IF(NC.GT.O)  GO  TO  13 
WRITE(6, 1 10) 

110  FORMAT  (1H0, ’ERROR  IN  DATA  INPUT') 

GO  TO  1000 

13  CONTINUE 

NDIM=0 
N=5+NC 
DO  1 4  1=1 , N 

IF(S(I) .NE.O.)  NDIM=NDIM+1 

14  CONTINUE 

DO  15  1=1, NC 

RW(I )=FLOAT(ND(NC+1 ) )/(FLOAT(ND(I+1 )-ND(I ) )*FLOAT(NC) ) 


ooooo  o  <r>  o  ooo 


C 


RW( I )= 1 . 0 
15  CONTINUE 


DO  16  1=1, N 
Y(I )=X(I ) 

16  CONTINUE 
C 

EST=1 . E-6 
EPS= 1 . E-6 
LIMIT=1 
IER=0 
C 


C***«*»*»***«*  FMCG  SEARCH  START 

C 

CALL  FMCG ( FUNCT , N , X , V , G , EST , EPS , LIMIT , IER , H ) 
««»«»«««  fMcg  SEARCH  •••••••••••••••*••**•**  END 


CALL  TITLE(N,Y, LIMIT, EPS) 

CALL  PRINT1(NfX,V,GtIER) 

CALL  NMBC(X , L , NC , WN , CF , PO , TO ) 

1000  CONTINUE 
STOP 
END 

SUBROUTINE  FUNCT(N , X , V , G) 

COMPUTATION  OF  THE  FUNCTION  VALUE  AND  DERIVATIVES 

(DOUBLE  EXPONENTIAL  FUNCTION) 

DIMENSION  X(9),6(9),F(9) 

COMMON  /PARM 1 /  NC , ND (5 ) , RW( 4 ) 

COMMON  /PARM2/  IC , PLOG (900 ) , TLOG (900 ) , ULOG(90O) , TAU(900 ) , S(9 ) 
C 

IC=IC+1 

V=0. 

DO  20  K=1,N 
G (K  )  =  0 . 

20  CONTINUE 
C 

DO  21  J  =  1  ,  NC 
JJ=J+5 
SQER=0 . 

DO  22  L  =  1 ,5 
F(L)=0. 

22  CONTINUE 
F(JJ)=0. 

IM=ND( J )+1 
IN=ND( J+1 ) 

DO  23  I =IM , IN 

W1 =X( J J )+X(4 )*PLOG(I )+X(5 )#TLOG(I )+ULOG(I ) 

R=X(1)+X(2)«WUX(3)*W1«W1 

R=10.*«R 


o  o  o 
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IFCR.LE.70.)  GO  TO  24 
TC=0. 

GO  TO  25 

24  CONTINUE 
TC=EXP(-R) 

25  CONTINUE 
E=TAU(I)-TC 
R=R*E*TC 
SQER=SQER+E**2 
F(  1  )  =F( 1 )+R 

F( 2 ) =F ( 2 )+R*W1 

F( 3 )=F( 3 )+R*W1 *W1  .  . 

R=R*(X(2)+2.*X(3)*W1) 

F(4)=F(4 )+R*PLOG ( I ) 

F(5)=F(5)+R*TLOG(I) 

F(JJ)=F(JJ)+R 
23  CONTINUE 

V=V+SQER*RW(  J) 

DO  26  K  =  1 ,5 
G(K)=G(K)+F(K)*RW(J) 

26  CONTINUE 
G(JJ)=F(JJ)*RW(J) 

21  CONTINUE 

DO  27  1  =  1  ,N 

G(I)=4.60517*G(I)*S(I) 

27  CONTINUE 

RETURN 

END 

SUBROUTINE  TITLE (N , X , LIMIT , EPS ) 

PRINTING  OF  THE  TITLE  AND  INITIAL  VALUES 

DIMENSION  X(9),L(4) 

COMMON  /PARM 1 /  NC , ND (5 )  ,  RW( 4 ) 

COMMON  /PARM4/  PO , TO , NDIM , ID (5 , 9 ) 

DO  40  1=1 ,NC 
L(I)=ND(I+1)-ND(I) 

40  CONTINUE 

CALL  DATE  ( MONTH , IDAY , IYEAR ) 

WRITE (6, 1 1 1 ) MONTH, IDAY, IYEAR 
111  F0RMAT(1H1,T60,I4, •  /  ',12,  1  /  ’,12,/) 

WRITE (6 , 400 )  NC , NDIM 

400  FORMAT  (1H  , T 1 4 , ' *»»  SIMULTANEOUS  PARAMETER  EVALUATION  »«*',///, 

*  •  PARAMETERS  :  (  N  ,  M  ,  A 1  ,  A2  ,  A3  ,  C(I),  1  =  1, ’,12,’  )', 

«  8X ,  '  (  DIMENSION  =’,I3,'  )',///,’  DATA  :') 

WRITE (6, 401)  ( ( ID (K , J ) ,J=1 ,9) ,L(K) ,K=1 ,NC) 

401  F0RMAT(1H  +  ,T1 1  ,  '  (  '  ,9A4,  ’  )',5X,'//  OF  POINTS  =',I5,//> 

WRITE (6, 402)  ND(NC+1 ), PO, TO , LIMIT , EPS 

402  FORMAT  ( 1H  +  .T51  ,  'TOTAL  //  OF  DATA  =',I5,///, 

»  '  FUNCTION  :  TAU  (  W  )  =  EXP  (  -10  ••  (  A1  +  A2  *  W  +' , 

»  •  A3  #  W**2  +  A4  *  W**3  )  )’ ,//,T15, 'WHERE,  ’, 


o  o  o  o  o  o  o 
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C 


C 


•  '  W  =  LOG(C)  +  LOG (  U  •  (P/PO)»»N  *  (TO/T)*«M  ) 

•  ' ,  A4  =  0. ' ,///, 

*  '  CONSTANTS  :  PO  = • ,F8.2,7X, 'TO  =',F8.2,7X, 

*  'LIMIT  =’ ,I5,7X, 'EPS  = • , 1PE10. 1 , /) 


403 

404 

405 

406 


IF(NC.EQ.I)  WRITE (6 , 403 )  X ( 1 ) , X( 6 ) , X (2 ) , X( 3 ) , X(4 ) , X(5 ) 

IF(NC . EQ. 2)  WRITE (6 , 404 )  X ( 1 ) ,X(6 ) , RW( 1 ) ,X(2 ) ,X(7 > , RW(2 ) , X( 3 ) , 

X(4),X(5) 

IF (NC . EQ . 3 )  WRITE (6 , 405 )  X ( 1 ) , X(6 ) , RW( 1 ) , X(2 ) , X (7 ) , RW(2 ) ,X( 3 ) , 

X(8),RW(3)fX(4),X(5) 

IF(NC . EQ. 4)  WRITE (6 , 406 )  X ( 1 ) , X(6 ) , RW( 1 ) , X(2 ) ,X(7 ) , RW(2 ) ,X( 3 ) , 

X(8) ,RW(3) ,X(4) , X ( 9 ) ,RW(4) ,X(5) 

FORMAT  (1H0, 'INITIAL  VALUES  :*,T22,'A1  =  ’  , F 12. 7 , 9X , ’LOG( C 1 )  =  ', 
F12.7.//.T22, 'A2  =' .FI  2 . 7 , // , T22 , ' A3  = ' ,F12 . 7 , // ,T22 , 'N  =  ’  , 

F12 . 7 , // , T22 , ' M  =',F12.7,/) 

FORMAT  (1H0, 'INITIAL  VALUES  :',T22,’A1  =  '  ,  F 1 2 . 7 , 9X  ,  ' LOG (Cl)  =  ', 
F12.7,4X,'(  WEIGHT  =’,F8.4,’  )',//, T22 ,' A2  =',F12.7,9X, 

'L0GCC2)  =' .F12.7.4X, ' (  WEIGHT  =’,F8.4,'  )',//, T22 ,' A3  =’, 
F12.7,// ,T22, 'N  = ' ,F 12 . 7 , // , T22 , ' M  «',F12.7,/> 

FORMAT  (1H0, 'INITIAL  VALUES  :',T22,’A1  = ' ,  F 1 2 . 7 , 9X  , ' LOG (Cl)  =’, 
F12.7,4X,’(  WEIGHT  =',F8.4,’  )',//, T22 A2  =',F12.7,9Xf 
'LOG ( C2)  =' ,F12.7,4Xf  '(  WEIGHT  =',F8.4,'  ) ’ , // , T22 , ' A3  =  ’, 
F12.7.9X, ’L0G(C3)  * ' ,F12. 7, 4X, ' (  WEIGHT  =',F8.4,'  )',//, T22, 

•N  ='  ,F 1 2 . 7 , // , T22 , 'M  =',F12.7,/) 

FORMAT  (1H0, 'INITIAL  VALUES  :',T22,’A1  = ' , F 1 2 . 7 , 9X f 'LOG ( C 1 )  r', 
FI  2 . 7 , 4X ,  '  (  WEIGHT  =\F8.4,'  )  • ,//,T22, »A2  s',F12.7,9X, 

' LOG ( C2 )  =' ,F12.7,4X, ' (  WEIGHT  =»,F8.4,»  ) • ,//,T22, 'A3  =’, 
F12.7.9X, 'L0GCC3)  = '  ,F 1 2 . 7 , 4X , ' (  WEIGHT  =',F8.4,’  )’,//, T22, 

'N  =  '  ,  F 12 . 7 , 9X , 'LOG (C4)  = ’  ,F 12 . 7 , 4X , ' (  WEIGHT  =',F8.4,'  )’,//, 
T22,'M  =  '  , F 12 . 7 , /) 


RETURN 

END 

SUBROUTINE  PRINT  1 (N , X , V , G , IER ) 


PRINTING  OF  THE  RESULTS  AND  COMPUTATION/PRINTING  OF  ERRORS 
AND  STANDARD  DEVIATIONS 


DIMENSION  X(9) ,G(9) ,E(900) ,PD(900) ,TC(900) ,W(900) 

COMMON  /PARM 1 /  NC , ND(5 )  ,  RW(4 ) 

COMMON  /PARM2/  IC , PLOG ( 900 ) , TLOG (900 ) , ULOG ( 900 ) , TAU(900 ) , S(9 ) 
COMMON  /PARM3/  P (900 ) , T (900 ) , U ( 900 ) , L (20 ) 

COMMON  /PARM4/  PO , TO , NDIM , ID (5 , 9  ) 

EQUIVALENCE  (E,PLOG) , (PD, TLOG) , (TC.ULOG) 

EQUIVALENCE  IS  FOR  SPACE  CONSERVATION 


ITsO 

TTD=0. 

TV=0 . 

V=SQRT (V /FLOAT (ND(NC+1))) 

WRITEC6 ,510)  IER,IC,(X(I),G(I),I=1,5) 

510  FORMAT  (1H0,'»*  RESULTS  OF  COMPUTATION  :  IER  =',I3,4X, 

•  'SUBROUTINE  FUNCT  WAS  CALLED', 16,'  TIMES  «»', 

•  ///,'  FINAL  VALUES  AND  GRADIENTS  :', 
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C 


C 


c 


c 


*  T35, ' A 1 

*  T35, '  A2 

*  T35,'A3 

*  T35 ,  ' N 

*  T35,'M 
WRITE(6 ,511) 


=',F12.7,T65,'D/D(A1) 

=  '  ,F12.7,T65, 'D/D(A2) 

=  '  ,F12.7,T65, 'D/D(A3) 

=  ' ,F  12 . 7 , T65 , ' D/D (N ) 

=  ' ,F12.7,T65, ’D/D(M) 

(I  ,X(I+5) , I , G(I+5 ) ,1  =  1 ,NC) 


=  ' , E15. 6, // , 
=  ’ ,E15.6,//, 
=' ,E15.6,//f 
=  V, E15 . 6 , // , 
=  ’ , E15 . 6 ) 


511  FORMAT ( 1  HO, T35, 'LOG(C' , 1 1 , • )  =  * ,F 12 . 7 , T65 , 'D/D(LOG(C' ,11 ,'))=* , 


*  E15.6) 

WRITE(6 ,512)  V 

512  FORMAT(/,1HO, 'FINAL  STANDARD  DEVIATION  : * ,F 1 5 . 7  ) 


WRITE (6 ,513)  (L ( I )  ,1=1 ,20) 

513  FORMAT(1HO,T5,  ’  (  COMMENT  :  \20A4,'  )') 


DO  50  J=1 ,NC 
J J=J+5 
V  =  0. 

TD=0 . 

IM=ND( J )+1 
IN=ND( J+1 ) 

K=ND( J+1 )-ND( J ) 
RK=FLOAT (K ) 


DO  51  I =IM , IN 

W(I)=X(JJ)+X(4) *PLOG( I )+X (5 )*TLOG( I )+ULOG( I ) 
R=10.**(X(1)+X(2)*W(I)+X(3)*W(I)**2) 

IF ( R . LE . 70 )  GO  TO  52 
TC ( I ) =0 . 

GO  TO  53 

52  CONTINUE 
TC(I)=EXP(-R) 

53  CONTINUE 
E(I)=TAU(I)-TC(I) 

PD(I)=100.*E(I)/TAU(I) 

TD=TD+ABS(E ( I ) ) 

V=V+E ( I ) **2 

51  CONTINUE 

TTD=TTD+TD 
TV=TV+V 
TD=TD/RK 
V=V/RK 
SF =SQRT ( V ) 

WRITE (6 ,514)  (ID(J,I) ,1=1 ,9),K, (X ( I ) , 1= 1 , 5 ) , J , X( J J ) , TD, V, SF 

514  FORMAT  ( 1 H 1 , T 1 5 ,' ACTUAL/COMPUTED  TRANSMITTANCES ',//// , 

*  T 10 , ' DATA  : ' , 9A4, 5X, '#  OF  POINTS  = ' , 14, //// , T10 , 

*  ' A  1  =' ,F12.7,8X, 'A2  = • ,F12 .7, 8X, 'A3  = ' ,F 12 . 7 , // , T1 0 , 'N  =',F12.7, 

»8X,'M  =' ,F12.7,8X, 'LOG(C*  ,11 ,  '  )  = ' ,F 1 2 . 7 , //// , 

*  T15, 'AVERAGE  DISCREPANCY  = ' , F 12 . 7 , // , 

*  T15 , 'MEAN  SQUARE  ERROR  =',F12.7,//, 

*  T15, 'STANDARD  DEVIATION  = ' ,F 1 2 . 7 , //// , 

«  T5, '#' ,T11, 'U  = ' , T26 , ' P  =',T36,'T  =',T46,'X  =', T56 ,' ACTUAL ’ , 

»  T65, 'COMPUTED' ,T76, 'DIFFERENCE' ,T90, '*  DIFF.',/) 

WRITE(6 ,515)  (I,U(I),P(I),T(I),W(I),TAU(I),TC(I),E(I),PD(I), 

«  I  =  IM , IN  ) 

515  FORMAT  (1H  , 15 , T10 , El  1 . 4 , T24 ,F8 . 2 , T34 ,F8 . 2, T43 ,F9 . 4 , T53 ,F9 . 4 , T63, 
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*  F9.4,T74,F11 .6,T87,F10.4) 

50  CONTINUE 

FK=FL0AT(ND(NC+1 ) ) 

TTD=TTD/FK  .  . 

TV=TV/FK 
TSD=SQRT (TV) 

WRITE(6 ,516)  NDCNC+1 ) ,TTD,TV,TSD 

516  FORMAT  ( 1H 1 ,///, T10 , 'TOTAL  If  OF  POINTS  USED  =  •  ,  16 , /// ,  T1 0  , 

*  'GLOBAL  AVERAGE’  DISCREPANCY  ='  , F 1 2 . 7 , /// , T1 0 , 

«  ’GLOBAL  MEAN  SQUARE  ERROR  = ’ , F 12 . 7 , /// , T 10 , 

*  'GLOBAL  STANDARD  DEVIATION  =’,F12.7) 

RETURN 

END 

SUBROUTINE  NMBC ( X , NAME , NC , WN , CF , PO , TO ) 

COMPUTATION  OF  THE  C '-VALUES  FOR  NON-MAJOR  BANDS 

DIMENSION  X(  9  )'  i  NAME (20  )  ,WN  ( 4  ) 

DIMENSION  CS(10) ,FS(10) 

DF=1 .E30 
SGN  =  1  . 

IF(X(3).LT.O.)  SGNs-1 . 

IF  THE  QUADRATIC  TERM  IS  TOO  SMALL,  THEN  IT  WILL  BE  IGNORED 

SMI=-2 . *X ( 3 ) /X (2 ) 

IF ( ABS (SMI ) .LE . 1 . E-6 )  GO  TO  50 
SYM=1 ./SMI 
50  CONTINUE 

WRITE(6,5)(NAME(I) ,1=1 ,20) 

5  FORMAT ( 1H1 ,T15,20A4) 

WRITEC6, 10) 

10  FORMAT ( 1H0, T15 , '  ***  CALCULATION  OF  THE  SPECTRAL  PARAMETER  FOR', 
«  '  NON-MAJOR  BANDS  *»*’,///) 

11  CONTINUE 
NFREQ=0 

12  CONTINUE  ’ 

C=0. 

1=0 

15  CONTINUE 

READ (5 , 20 , END=40 )  KGAS , FREQ , P , T, UGAS , TX 
20  FORMAT( 12, F 1 0 . 3 , E 1 1 . 4 ,F9 . 3 , 24X , E 1 1 . 4 ,F7 . 4 ) 

IF (KGAS. EQ. 0)  GO  TO  25 
IF(KGAS.LT.O)  GO  TO  35 
IF(UGAS.GE.DF)  GO  TO  15 


1  =  1  +  1 

WX=FREQ 

UGAS=UGAS/CF 

IF (SMI . LE . 1 . E-6 )  GO  TO  51 

CASE  1  QUADRATIC  TERM  IS  LARGE  AND  USED 


Mm 


o  o  o 


180 


XS=SYM+SGN*ABS (SQRT (X(2)**2-4.*X(3)*(X(1 ) -A LOG  1 0 (-ALOG (TX ) ) ) ) 
*  /( 2 .  *X ( 3 )  ) ) 

GO  TO  52 

51  CONTINUE 

CASE  2  QUADRATIC  TERM  IS  SMALL  AND  IGNORED 

XS= ( ALOG 10 (-ALOG (TX) )-X(1))/X(2) 

52  CONTINUE 

XC=X (4 )* ALOG 10 (P/PO )+X (5 ) *ALOG 10 (T0/T)+ALOG10(UGAS) 
C=C+(XS-XC) 

GO  TO  15 

25  C=C/FLOAT ( I ) 

NFREQ=NFREQ+1 
CS (NFREQ) =C 
FS(NFREQ)=WX 
DO  27  M=1,NC 

IF ( ABS ( WX-WN ( M) ) .LE.O. 1 )  CS (NFREQ) =X (5+M) 

27  CONTINUE 

IF (NFREQ. EQ. 10)  GO  TO  30 
GO  TO  12 

30  CONTINUE 

WRITE (6, 31 )(FS(K) ,K=1 , NFREQ) 

31  F0RMAT(1H0,2X, ’WAVE  NUMBER 2X , 1  OF  1 1 . 0 ) 

WRITE (6, 32 )( CS (K ) ,K=1 .NFREQ) 

32  FORMAT ( 1H0 , 5X , ’C  VALUES ' ,2X, 10F1 1 . 3//> 

GO  TO  11 

35  CONTINUE 

IF (NFREQ. EQ.  0)  GO  TO  40 
WRITE (6, 31 )(FS(K) ,K=1 .NFREQ) 

WRITE (6, 32 )( CS (K ) ,K=1 .NFREQ) 

40  CONTINUE 
RETURN 
END 


*1 


